{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Secure Kaplan-Meier Survival Analysis Explained\n",
    "\n",
    "The MPyC demo [kmsurival.py](kmsurvival.py) implements privacy-preserving Kaplan-Meier [survival analysis](https://en.wikipedia.org/wiki/Survival_analysis), based on earlier work by Meilof Veeningen. The demo is built using the Python package [lifelines](https://pypi.org/project/lifelines/), which provides extensive support for survival\n",
    "analysis and includes several datasets. We use lifelines for plotting Kaplan-Meier survival curves, performing logrank tests to compare survival curves, and printing survival tables."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "# use pip (or, conda) to make sure lifelines package is installed:\n",
    "!pip -q install lifelines  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "import os, functools\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import lifelines.statistics\n",
    "from mpyc.runtime import mpc\n",
    "mpc.logging(False)\n",
    "from kmsurvival import fit_plot, events_to_table, events_from_table, logrank_test, aggregate, agg_logrank_test"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Actual use of the lifelines package is hidden mostly inside the [kmsurival.py](kmsurvival.py) demo, except for the function `lifelines.statistics.logrank_test()` used below to validate the results of the secure computations.\n",
    "\n",
    "## Kaplan-Meier Analysis\n",
    "\n",
    "We analyze the aml (\"acute myelogenous leukemia\") dataset, which is also used as a [small example on Wikipedia](https://en.wikipedia.org/wiki/Survival_analysis#Example:_Acute_myelogenous_leukemia_survival_data). The file `aml.csv` (copied from the R package `KMsurv`) contains the raw data for 23 patients. Status 1 stands for the event \"recurrence of aml cancer\" and status 0 means no event (\"censored\")."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<style  type=\"text/css\" >\n",
       "</style><table id=\"T_2e9366d9_1b86_11eb_8735_989096b90405\" ><thead>    <tr>        <th class=\"col_heading level0 col0\" >observation</th>        <th class=\"col_heading level0 col1\" >time</th>        <th class=\"col_heading level0 col2\" >status</th>        <th class=\"col_heading level0 col3\" >group</th>    </tr></thead><tbody>\n",
       "                <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row0_col0\" class=\"data row0 col0\" >12</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row0_col1\" class=\"data row0 col1\" >5</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row0_col2\" class=\"data row0 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row0_col3\" class=\"data row0 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row1_col0\" class=\"data row1 col0\" >13</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row1_col1\" class=\"data row1 col1\" >5</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row1_col2\" class=\"data row1 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row1_col3\" class=\"data row1 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row2_col0\" class=\"data row2 col0\" >14</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row2_col1\" class=\"data row2 col1\" >8</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row2_col2\" class=\"data row2 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row2_col3\" class=\"data row2 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row3_col0\" class=\"data row3 col0\" >15</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row3_col1\" class=\"data row3 col1\" >8</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row3_col2\" class=\"data row3 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row3_col3\" class=\"data row3 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row4_col0\" class=\"data row4 col0\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row4_col1\" class=\"data row4 col1\" >9</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row4_col2\" class=\"data row4 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row4_col3\" class=\"data row4 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row5_col0\" class=\"data row5 col0\" >16</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row5_col1\" class=\"data row5 col1\" >12</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row5_col2\" class=\"data row5 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row5_col3\" class=\"data row5 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row6_col0\" class=\"data row6 col0\" >2</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row6_col1\" class=\"data row6 col1\" >13</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row6_col2\" class=\"data row6 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row6_col3\" class=\"data row6 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row7_col0\" class=\"data row7 col0\" >3</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row7_col1\" class=\"data row7 col1\" >13</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row7_col2\" class=\"data row7 col2\" >0</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row7_col3\" class=\"data row7 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row8_col0\" class=\"data row8 col0\" >17</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row8_col1\" class=\"data row8 col1\" >16</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row8_col2\" class=\"data row8 col2\" >0</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row8_col3\" class=\"data row8 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row9_col0\" class=\"data row9 col0\" >4</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row9_col1\" class=\"data row9 col1\" >18</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row9_col2\" class=\"data row9 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row9_col3\" class=\"data row9 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row10_col0\" class=\"data row10 col0\" >5</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row10_col1\" class=\"data row10 col1\" >23</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row10_col2\" class=\"data row10 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row10_col3\" class=\"data row10 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row11_col0\" class=\"data row11 col0\" >18</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row11_col1\" class=\"data row11 col1\" >23</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row11_col2\" class=\"data row11 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row11_col3\" class=\"data row11 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row12_col0\" class=\"data row12 col0\" >19</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row12_col1\" class=\"data row12 col1\" >27</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row12_col2\" class=\"data row12 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row12_col3\" class=\"data row12 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row13_col0\" class=\"data row13 col0\" >6</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row13_col1\" class=\"data row13 col1\" >28</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row13_col2\" class=\"data row13 col2\" >0</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row13_col3\" class=\"data row13 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row14_col0\" class=\"data row14 col0\" >20</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row14_col1\" class=\"data row14 col1\" >30</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row14_col2\" class=\"data row14 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row14_col3\" class=\"data row14 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row15_col0\" class=\"data row15 col0\" >7</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row15_col1\" class=\"data row15 col1\" >31</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row15_col2\" class=\"data row15 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row15_col3\" class=\"data row15 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row16_col0\" class=\"data row16 col0\" >21</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row16_col1\" class=\"data row16 col1\" >33</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row16_col2\" class=\"data row16 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row16_col3\" class=\"data row16 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row17_col0\" class=\"data row17 col0\" >8</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row17_col1\" class=\"data row17 col1\" >34</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row17_col2\" class=\"data row17 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row17_col3\" class=\"data row17 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row18_col0\" class=\"data row18 col0\" >22</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row18_col1\" class=\"data row18 col1\" >43</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row18_col2\" class=\"data row18 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row18_col3\" class=\"data row18 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row19_col0\" class=\"data row19 col0\" >9</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row19_col1\" class=\"data row19 col1\" >45</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row19_col2\" class=\"data row19 col2\" >0</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row19_col3\" class=\"data row19 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row20_col0\" class=\"data row20 col0\" >23</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row20_col1\" class=\"data row20 col1\" >45</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row20_col2\" class=\"data row20 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row20_col3\" class=\"data row20 col3\" >2</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row21_col0\" class=\"data row21 col0\" >10</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row21_col1\" class=\"data row21 col1\" >48</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row21_col2\" class=\"data row21 col2\" >1</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row21_col3\" class=\"data row21 col3\" >1</td>\n",
       "            </tr>\n",
       "            <tr>\n",
       "                                <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row22_col0\" class=\"data row22 col0\" >11</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row22_col1\" class=\"data row22 col1\" >161</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row22_col2\" class=\"data row22 col2\" >0</td>\n",
       "                        <td id=\"T_2e9366d9_1b86_11eb_8735_989096b90405row22_col3\" class=\"data row22 col3\" >1</td>\n",
       "            </tr>\n",
       "    </tbody></table>"
      ],
      "text/plain": [
       "<pandas.io.formats.style.Styler at 0x26f0a32dd00>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "df = pd.read_csv(os.path.join('data', 'surv', 'aml.csv')).rename(columns={'Unnamed: 0': 'observation', 'cens': 'status'})\n",
    "df.sort_values(['time', 'observation']).style.hide_index()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Time is in weeks. The study compares the time until recurrence among two groups of patients. Patients in group 1 received maintenance chemotherapy, while patients in group 2 did not get any maintenance treatment. To plot the [Kaplan–Meier survival curve](https://en.wikipedia.org/wiki/Kaplan–Meier_estimator) for groups 1 and 2, we use the function `fit_plot()` as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAakAAAGcCAYAAABjm19FAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAABFQklEQVR4nO3deZhU1bX+8e/LIKggKogaQEHjhCKgKCqKEByj0URNFM2NJDH+zE2ccvVepzglXpOYSW9MjEMkGpynEMVZiAaHgIoDoAYRFVREIpOCAq7fH2c3Fk01Xd10V52m3s/z1EPVGfZZdZqu1fucXWsrIjAzM8ujVpUOwMzMrC5OUmZmlltOUmZmlltOUmZmlltOUmZmlltOUmZmlltOUlY1JA2RNLPScayOpPslHV/pOBpK0iJJWzVBOyHpi00Rk60dnKQsdyTNkLRfwetjJH0oad9KxrU6kkakD9jf1Fp+eFo+spR2IuLgiPhzswTZjCKiQ0RMr3QctvZxkrJcS72KK4FDIuLvlY6nHq8D35DUpmDZ8cBrzX3gWsdsMW1XgqTWlY7BSuckZbkl6f8BvwIOjIgn07KtJT0maa6kDySNkrRhwT4zJJ0taUrqfV0vqX0d7Z8l6XVJC9P2XytYN0LSPyT9MrXzhqSD6wn5PeAl4MDUxsbAXsDoWsfdQ9KTkuZJekHSkIJ14ySdUPD6O5KmphgelLRlwbqQ9ANJ/wL+VeT9tZf0l3Su5kmaIGnTgvNU2Fu9UNJf0vOeqe3vSnoLeCxdhvxhrfZfkHREQSxflDRQ0nuFiUDS1yS9mJ7vLumpFM+7kn4naZ16zmtNOxunn+c76Xzck5aPkPSPWtuuuGwoaaSkP0gaI+kj4Ix6YmxV8H9jrqTb0s9ytefUmoeTlOXV94GLgWERMbFguYBLgS8AOwA9gAtr7XscWaLYGtgWOK+OY7wO7AN0Ai4C/iJp84L1A4FXgS7AL4DrJKmeuG8AvpWeHwP8FfhkRfBSN+A+4KfAxsAZwJ2SNqndkKTDgXOAI4BNgCeAm2tt9tUUZ+8isRyf3lsPoDNwErC4nvgL7Ut2jg9Mxx1eEFtvYMv0XlaIiGeAj4AvFSw+FrgpPV8OnE52TvcEhgH/WWI8NwLrATsCXYHfrH7zlRwLXAJ0BC6vJ8aTyc7rvmT/zz4k683Dmp9TayAnKcur/YGnyXomK0TEtIh4OCI+iYg5wK/JPkwK/S4i3o6If5N9MA2niIi4PSLeiYjPIuJWst7I7gWbvBkR10TEcuDPwOZAfX813w0MkdSJLFndUGv9N4ExETEmHfdhYCLw5SJtnQRcGhFTI2IZ8L9Av8LeVFr/74go9kG5lOyD9IsRsTwino2IBfXEX+jCiPgotX13rWMfB9wVEZ8U2W9FQpPUMb23mwFSDE9HxLKImAH8kVV/fqtIfzwcDJwUER9GxNIGXv79a0SMT+d8yepiJDvv50bEzPT+LgSOSpc91/ScWgM5SVlefZ+sF3RtYe9F0qaSbpE0S9IC4C9kf5UXervg+Ztkfw2vQtK3JE1Kl23mATvVauu9micR8XF62kHSPspGsy2SNLmwzfSBfh9Z761zRIyvddgtga/XHDMdd2+yBFjblsDlBdv9m6wn2a2O91rbjcCDwC3pEtkvJLVdzfa1rWg7Iham93VMWjQcGFXHfjcBR0hqR9YLfC4i3gSQtK2ke9PltgVkibf2z6+YHsC/I+LDBsRfqPZ5qjNGsvN+d8F5n0rWA9yUNT+n1kBOUpZXs8kuBe0D/L5g+f8CAfSJiA3Ieia1L8H1KHi+BfBO7cZTj+Aa4IdkyWRD4OUiba0iIp5Io9k6RMSORTa5AfgvsgRa29vAjRGxYcFj/Yj4WR3b/r9a265bc3+uJpzVxLk0Ii6KiN5k98YO5fNLkR+RXTqrsVmxJmq9vhkYLmlPoD0wto7jTiH74+BgVr6MBvAH4BVgm/TzO4cSzjnZudhYBfcfC6z0XiTV+17qifFt4OBa5719RMyq55xaM3CSstyKiHfIEtVB+nxod0dgETA/3d85s8iuP5DUPd3sPhe4tcg265N9cM0BkPRtsp5UU/g72eXK/yuy7i/AVyQdKKl1uhE/RFL3ItteBZwtaccUYydJXy81CElDJfVJAwQWkF2q+iytngQcI6mtpAHAUSU0OYasl3ExcGtEfLaabW8CTgUGA7cXLO+YYlkkaXuyHnO9IuJd4H7g95I2SnEPTqtfAHaU1E/ZIJkLS2lzNTFeBVxSc2lT0ibp/mB959SagZOU5VpEvEV2g/soSZeSDXDYBZhPdvnpriK73QQ8BEwnGxzx0yLtTiEbOfgUWa+tD1D70lxjY46IeDTdE6u97m2gZkDEHLK/2s+kyO9iRNwN/Jzs0tICsp5efSMMC20G3EH2YTqVLHnemNb9mGxgyYdk5/SmYg3UiucTsvO9Xwnb30x2r+mxiPigYPkZZD2XhWQ92WJ/QNTlP8iSwivA+8BpKa7XyBLnI2T3Ff9Rx/6lxng52YjMhyQtJLs3OjCtW905tWYgT3poaxNJM4ATIuKRSsdiZmvOPSkzM8stJykzM8stX+4zM7Pcck/KzMxyy0nKzMxya62qbpwHXbp0iZ49e1Y6DDOzFuXZZ5/9ICJWqWHpJNXEevbsycSJE+vf0MzMVpD0ZrHlvtxnZma55SRlZma55SRlZma55XtSZlZ2S5cuZebMmSxZsqTSoViZtW/fnu7du9O2bWkznFRtkpL0J7Iy++9HxCrVr9McRpeTTYb2MTAiIp4rb5Rma6eZM2fSsWNHevbsSf2THdvaIiKYO3cuM2fOpFevXiXtU82X+0YCB61m/cHANulxItk8OGbWBJYsWULnzp2doKqMJDp37tygHnTVJqmIeJxsptO6HA7ckKZdeBrYME1hbWZNwAmqOjX05161l/tK0I2Vp5yemZa92xwHe/r336PjvKklbTt+3aE8ut6Xi647vF83jh24RVOGZrZW+s53vsO9995L165defnll0vaZ9y4cQwdOpRrrrmGE044AYBJkybRv39/LrvsMs4444w6973qqqtYb731+Na36p7Id9KkSbzzzjt8+cvFf79rTJw4kRtuuIErrriipLhXZ+TIkUycOJHf/e53a9xWc6janlRTknSipImSJs6ZM6fR7Xy2bFm9jy0/fZ09Fz7MgiVLV3m8NGs+Nz49o+nemNlabMSIETzwwAMN3m+nnXbitttuW/H65ptvpm/fvvXud9JJJ602QUGWpMaMGVNvWwMGDGiSBNUSuCdVt1lAj4LX3dOyVUTE1cDVAAMGDGhUWfk9/vOa0ja8/hD6LJnH/d8fvMqqo//4FAuWLG3M4c2qzuDBg5kxY0aD99tyyy1ZsGABs2fPpmvXrjzwwAMr9XyuueYarr76aj799FO++MUvcuONN7Leeutx4YUX0qFDB8444wyGDBnCwIEDGTt2LPPmzeO6665j4MCBnH/++SxevJh//OMfnH322fTq1YtTTz2VJUuWsO6663L99dez3XbbMW7cOH75y19y7733cuGFF/LWW28xffp03nrrLU477TROOeUUAP7yl79wxRVX8OmnnzJw4EB+//vf07p1a66//nouvfRSNtxwQ/r27Uu7du2a6rQ2OSepuo0GfijpFrKpo+dHRLNc6mtKyz8LRk8qmksB6LhuW4Zu17WMEZmt3kV/m8yUdxY0aZu9v7ABF3xlxwbvd9lllzFq1KhVlg8ePHilnstRRx3F7bffTv/+/dlll11W+pA/4ogj+N73vgfAeeedx3XXXcfJJ5+8SpvLli3jn//8J2PGjOGiiy7ikUce4eKLL17p0tuCBQt44oknaNOmDY888gjnnHMOd9555yptvfLKK4wdO5aFCxey3Xbb8f3vf59p06Zx6623Mn78eNq2bct//ud/MmrUKPbff38uuOACnn32WTp16sTQoUPp379/g89VuVRtkpJ0MzAE6CJpJnAB0BYgIq4CxpANP59GNgT925WJtGE++yzYpGP7OtfPWejvpZjV5cwzz+TMM8+sd7tvfOMbHH300bzyyisMHz6cJ598csW6l19+mfPOO4958+axaNEiDjzwwKJtHHHEEQDsuuuudfbo5s+fz/HHH8+//vUvJLF0afErJYcccgjt2rWjXbt2dO3aldmzZ/Poo4/y7LPPsttuuwGwePFiunbtyjPPPMOQIUPYZJOsluvRRx/Na6+9Vu97rpSqTVIRMbye9QH8oEzhmFWtxvR4mkupPanNNtuMtm3b8vDDD3P55ZevlKRGjBjBPffcQ9++fRk5ciTjxo0reqya3lfr1q1ZtmxZ0W1+/OMfM3ToUO6++25mzJjBkCFDVttWYXsRwfHHH8+ll1660rb33HNP0TbyqmqTlJlZbaX2pAAuvvhi3n//fVq3br3S8oULF7L55puzdOlSRo0aRbdu3Uo+fseOHVm4cOGK1/Pnz1+x/8iRI0tuB2DYsGEcfvjhnH766XTt2pV///vfLFy4kIEDB3Lqqacyd+5cNthgA26//faSBn5Uikf3mVlVGj58OHvuuSevvvoq3bt357rrrmvQ/nvttRdf/epXV1n+k5/8hIEDBzJo0CC23377BrU5dOhQpkyZQr9+/bj11lv57//+b84++2z69+9fZ2+rLr179+anP/0pBxxwADvvvDP7778/7777LptvvjkXXnghe+65J4MGDWKHHXZoULvlpuyqljWVAQMGRLPOJ3X9IbBkHnx//Cqrjv7jU8xd9Ak/+WqfOnefs3AJh/Ur/S87s+YwderU3H84WvMp9vOX9GxEDKi9rS/3tUSxHF66Y+Vl7TsBHSoSjplZc3GSaok++ww6bLryskWzKSVJLVm2fLVD1MHD1M0sP5ykqkyPjdavdxsPUzezvPDACTMzyy0nKTMzyy0nKTMzyy0nKTOrOm+//TZDhw6ld+/e7Ljjjlx++eUl7ztu3Dgk8be//W3FskMPPbTOyhI1fvvb3/Lxxx83NuTV2muvverdptTjn3DCCUyZMqUpwqJnz5588MEHa9SGk5SZVZ02bdrwq1/9iilTpvD0009z5ZVXNuiDuXv37lxyySUNOmZzJqnCskxrevxrr72W3r17N0VYTcJJysyqzuabb84uu+wCZKWIdthhB2bNWv1XMwr17duXTp068fDDD6+y7tFHH6V///706dOH73znO3zyySdcccUVvPPOOwwdOpShQ4eusk/Pnj05++yz6devHwMGDOC5557jwAMPZOutt+aqq64CYNGiRQwbNoxddtmFPn368Ne//nXF/h06ZF8/GTduHEOGDOGoo45i++2357jjjiMiih7/+9//PgMGDGDHHXfkggsuWNHWkCFDqClI0KFDB84991z69u3LHnvswezZswGYM2cORx55JLvtthu77bYb48dnxQXmzp3LAQccwI477sgJJ5xAUxSL8BD0tcyseYu5+N7Jda4ftHUXhu2waZ3rzcru/rPgvZeats3N+sDBPytp0xkzZvD8888zcOBAoPQis+eeey4//vGP2X///VcsW7JkCSNGjODRRx9l22235Vvf+hZ/+MMfOO200/j1r3/N2LFj6dKlS9E4tthiCyZNmsTpp5/OiBEjGD9+PEuWLGGnnXbipJNOon379tx9991ssMEGfPDBB+yxxx4cdthhq0zH/vzzzzN58mS+8IUvMGjQIMaPH88pp5yyyvEvueQSNt54Y5YvX86wYcN48cUX2XnnnVdq66OPPmKPPfbgkksu4b//+7+55pprOO+88zj11FM5/fTT2XvvvXnrrbc48MADmTp1KhdddBF77703559/Pvfdd1+DS00V4yS1Fjm8XzfmLvqkzvVvzv0Y+MBJyixZtGgRRx55JL/97W/ZYIMNgNKLzA4enE08+o9//GPFsldffZVevXqx7bbbAnD88cdz5ZVXctppp9Xb3mGHHQZAnz59WLRoER07dqRjx460a9eOefPmsf7663POOefw+OOP06pVK2bNmsXs2bPZbLPNVmpn9913p3v37gD069ePGTNmsPfee69yvNtuu42rr76aZcuW8e677zJlypRVktQ666zDoYceCmRTitT0HB955JGVLo8uWLCARYsW8fjjj3PXXXcB2fQhG220Ub3vuz5OUmuRYwduQYd2reucT2p1PSyziimxx9PUli5dypFHHslxxx23Ym4nKL0nBVlv6qc//Slt2qz5R2nNdButWrVaaeqNVq1asWzZMkaNGsWcOXN49tlnadu2LT179mTJklW/eF9s2o7a3njjDX75y18yYcIENtpoI0aMGFG0rbZt267oqRW29dlnn/H000/Tvn3dc9c1FScpW0UppZPA5ZOs5YoIvvvd77LDDjvwox/9aKV1DZmu44ADDuDHP/4x776bTdq93XbbMWPGDKZNm7Zi6vh9990X+Hwajrou99Vn/vz5dO3albZt2zJ27FjefPPNBu1fePwFCxaw/vrr06lTJ2bPns39999f51xVxRxwwAH83//934rzNGnSJPr168fgwYO56aabOO+887j//vv58MMPGxRjMR44YavosdH6bNKxfb2PhYuLzxJqlnfjx4/nxhtv5LHHHqNfv37069ePMWPGNKqtc889l7fffhuA9u3bc/311/P1r3+dPn360KpVK0466SQATjzxRA466KCiAydKcdxxxzFx4kT69OnDDTfc0OBpQAqP37dvX/r378/222/Psccey6BBgxrU1hVXXMHEiRPZeeed6d2794rBHRdccAGPP/44O+64I3fddRdbbLFFg9otxlN1NLGyTNXx0Rw45FcrL180G/ocxehJs+q93Hf+oU0zE6qn/bDG8lQd1a0hU3W4J2VmZrnlJGVmZrnlJGVmZrnlJGVmFeH74dWpoT93D0FfWyxdAi/dwRfe/pAN1m27yupl62wANG7oq1lTa9++PXPnzqVz586rVEywtVdEMHfu3AZ9v8pJam2x0ZYAfNK+NUvXXWeV1W0Xz8FJyvKie/fuzJw5kzlz5lQ6FCuz9u3br6iIUQonKTMru7Zt29KrV69Kh2EtgO9JVZk5C+uu7WdmljdOUlXmg0WfVjoEM7OS+XJfSzT/LXjgrJWXbTUEtj2orGGUWuOvGNf9M7NSOEm1NH2OysoiFfr3G9m/ZU5SPTZav9H7zlm4asVlM7PanKRamgHfhnYdoUPBnFC1e1VmZmsJJ6kqVN+8Up6918zyoqqTlKSDgMuB1sC1EfGzWuu3AP4MbJi2OSsiGlfPv8xGvQo3vVa4ZBNgIQBT3124yvZdOqzDJh3befZeM8uVqk1SkloDVwL7AzOBCZJGR8SUgs3OA26LiD9I6g2MAXqWPdhGOG677FGj7eI5zO31FYZf8zQ3f2+POvfz7L1mlifVPAR9d2BaREyPiE+BW4DDa20TwAbpeSfgnTLGZ2ZW9aq2JwV0A94ueD0TGFhrmwuBhySdDKwP7FesIUknAicCTTITZXNotXwJnd/4G7BJ+re4tos7Aax2mxrL1tmA+d32baoQzcxWUc1JqhTDgZER8StJewI3StopIj4r3Cgirgauhmxm3grEWa9POvRY8XzpupvUuV20rn+bGlk9QDOz5lPNl/tmAT0KXndPywp9F7gNICKeAtrjKq1mZmVTzUlqArCNpF6S1gGOAUbX2uYtYBiApB3IklSL7j4cu22lIzAzK13VJqmIWAb8EHgQmEo2im+ypIslHZY2+y/ge5JeAG4GRkQLn6mtcMSfmVneVfU9qfSdpzG1lp1f8HwKMKjccZmZWaZqe1JmZpZ/TlJmZpZbTlJmZpZbVX1Pyiqn2FxUnmPKzGpzkrKKKDYXleeYMrPafLnPzMxyy0nKzMxyy0nKzMxyy/ekrNE+r6y+KldIN7Om4CRljVZYWb02V0g3s6bgJGWrmD4fznqy/u327QYHb9n88ZhZ9XKSspXs26207abPz/51kjKz5uQkZSs5eMvSEk8pPS0zszXl0X1mZpZbTlJmZpZbTlJmZpZbTlJmZpZbTlJmZpZbTlJmZpZbTlJmZpZb/p7U2uLfb8ADZ7HT4qW0ab3q3x7zN9uLD7sPq0BgZmaN5yS1NthqyGpXt1/4JoCTlJm1OE5Sa4NtD8oewMvTPmDD9dZZaXXPiT+pRFQNVmxK+ebiqerNWgYnKcuNYlPKNxdPVW/WMjhJWbOoa64pzzNlZg3hJGXNomauqVGvwnHbfb7c80yZWUN4CLo1q5teq3QEZtaSOUmZmVluOUmZmVluOUmZmVluVfXACUkHAZcDrYFrI+JnRbb5BnAhEMALEXFsWYPMsenzS5uht3AbLe/E0smTGbR1F4btsGnzBWdma4WqTVKSWgNXAvsDM4EJkkZHxJSCbbYBzgYGRcSHktbKb39u8vodzNn6qAbts2+3lV/P/hjeX1x825fmFr5aB+YtZM7CT5ykzKxeVZukgN2BaRExHUDSLcDhwJSCbb4HXBkRHwJExPtlj7IMuk6/q8FJ6uAts0d9Dvkb3PeVz1+3XTyHUydv1cAIzaxaVfM9qW7A2wWvZ6ZlhbYFtpU0XtLT6fKgmZmVSTX3pErRBtgGGAJ0Bx6X1Cci5hVuJOlE4ESALbbYoswhmpmtvaq5JzUL6FHwuntaVmgmMDoilkbEG8BrZElrJRFxdUQMiIgBm2yySbMFbGZWbaq5JzUB2EZSL7LkdAxQe+TePcBw4HpJXcgu/00vZ5BNpf3CN1dbDb0xldI9R5WZNbeqTVIRsUzSD4EHyYag/ykiJku6GJgYEaPTugMkTQGWA2dGxNy6W82n+ZvtBWSDFtZZ8kHRbdb/cOoqyz5t34Wl6xbvGZY6R9Wx2zYkUjOzlVVtkgKIiDHAmFrLzi94HsCP0qPF+rD7sNUmkx0fPpbJ+9/UoDZL7XkVFpc1M2uoar4nZWZmOeckZWZmuVXVl/us/FotX0LbxdltvcJJET0ZopkV4yRlZfVJhx5E6+x54aAMT4ZoZsX4cp/x/lZHVDoEM7OinKSswXX7zMzKxUnKzMxyy0nKzMxyy0nKzMxyy0nKzMxyy0PQrSotWbac0ZNqF703szXRcd22DN2uaScwd5KyqtRjo/UrHYLZWmfOwiVN3qYv95mZWW65J7WWWa9dG+Z9/Gmd6z9d/hldO7YvY0RmZo3nJLWW6d9jw9WuHz+t+HxSZmZ55CRlFTF9Ppz15OevtbwTuy+ZzbAdNq1cUGaWO05S1mirm5J+dVPL79tt1WWvL2rD0tc/cJIys5U4SVmj1ExJX0x9U8sfvGX2KHT2E8tY2mTRmdnawkmqJWrfCRbNXnnZ0iWw0ZbFt28Gq5uSvtSp5c3M6uMk1RJts/+qy166o/xxmJk1M39PyszMcstJyszMcsuX+ywnPqPt4rl0fuNva9TKsnU2YH63fZsoJjOrNCcpy4VonVXBWLruJmvUTtvFc5oiHDPLCV/uMzOz3HJPqsrUV9sPXN/PzPLDSarK1FfbD1zfz8zyw5f7zMwst5ykzMwst5ykzMwst5ykzMwst6o6SUk6SNKrkqZJOms12x0pKSQNKGd8ZmbVrmqTlKTWwJXAwUBvYLik3kW26wicCjxT3gjNzKxqkxSwOzAtIqZHxKfALcDhRbb7CfBzYEk5gzMzs+pOUt2Atwtez0zLVpC0C9AjIu5bXUOSTpQ0UdLEOXNclsfMrKn4y7x1kNQK+DUwor5tI+Jq4GqAAQMGRPNG1jLUNbX86qaVbw53PPs2R+3ao2zHM7OmVc09qVlA4adX97SsRkdgJ2CcpBnAHsBoD56o3/zN9mJJx1VnCW6/8E06vfdkWWO587lZ9W9kZrlVzT2pCcA2knqRJadjgGNrVkbEfKBLzWtJ44AzImJimeNsceqaWt7TyptZQ1VtTyoilgE/BB4EpgK3RcRkSRdLOqyy0ZmZGVR3T4qIGAOMqbXs/Dq2HVKOmMzM7HNVnaQsX6bPh7PW8JaVlndi6eTJKy27+N7Jq2w3aOsuDNth0zU7mJk1Oycpy4V9u9W/zerM/hjeXwywDsxbuNK6qe8uXGX7OQs/cZIyawGcpCwXDt4ye6yptovnMLfXV1a8Hn7N09z8vT1W2qZYz8rM8qlqB06YmVn+uSdlqyhlinnI5zTzrZYvofMbfytYsknFYjGzNeckZasoZYp5yOc08590cHUJs7WJL/eZmVluuSe1tmjfCRbNrnv90iWwUROMTGhhvtnro0qHYGZrwElqbbHN/qtf/9Id5YkjZ7619cfMrXQQZtZovtxnZma55SRlZma55SRlZma55SRlZma55SRlZma55SRlZma55SHoVlbtF75ZdIbe+ZvtVXQ2XzOrbk5SVjbzN9ur6PL2C98EcJIys1U4SVnZfNh9WNFEVKxnZWYGvidlZmY55iRlZma55ct91aK+ArSNsN5H78J6WzVpm2ZmhZykqkV9BWgbYZ13/lTS5IjF5HHCRDPLHycpa7Qdv7ABdOjSqH3zOGGimeWP70mZmVluOUmZmVluOUmZmVluOUmZmVluOUmZmVluOUmZmVluOUmZmVluVXWSknSQpFclTZN0VpH1P5I0RdKLkh6VtGUl4jQzq1ZVm6QktQauBA4GegPDJfWutdnzwICI2Bm4A/hFeaM0M2s5xrz0bpO3WbVJCtgdmBYR0yPiU+AW4PDCDSJibER8nF4+DXQvc4xmZi3GA5Obtj4oVHeS6ga8XfB6ZlpWl+8C9zdrRGZmthLX7iuBpG8CA4B961h/InAiwBZbbFHGyNYedU0r31Ceht5s7VLNSWoW0KPgdfe0bCWS9gPOBfaNiE+KNRQRVwNXAwwYMCCaPtScWt30H0uXwEaljTOpa1r5BofjaejN1jrVnKQmANtI6kWWnI4Bji3cQFJ/4I/AQRHxfvlDzLnVTf/x0h0lN1PXtPIN1ZCe2JtzP+bieyev8THNrHlVbZKKiGWSfgg8CLQG/hQRkyVdDEyMiNHAZUAH4HZJAG9FxGEVC9qaxKCtuwCeKsSsseYs/IQPFhWfS67nWfetsuzUYdtw+v7bNupYVZukACJiDDCm1rLzC57vV/agrNkN22FThu2waaXDMFvrDL/maWb87JAmbbOaR/eZmVnOVXVPyipnvXZtVpl63lPKm1ltTlJWEf17bLjKMk8pb2a1+XKfmZnllpOUmZk1iYN2bPoBSU5SZmbWJL7cZ/Mmb9NJyszMcstJyszMcstJyszMcstD0K151FV8tgGFZ83MnKSsedRVfLYBhWfNzHy5z8zMcstJyszMcstJyszMcsv3pGytUnsa+kWdd2Zur69UMCIzWxNOUrbWqD0NffuFb6LlSysUjZk1BScpy41i03c0xLyN94GN91nxeqeXLqXVZ0uZs3DJKtsuWbacHhut3+hjmVl5OElZbhSbvmONTGsLy+Gwft1WWTV60qymPZaZNQsPnDAzs9xykjIzs9zy5T4rL5dLMrMGcJKy8nK5JDNrAF/uMzOz3HKSMjOz3HKSMjOz3HKSMjOz3HKSMjOz3HKSMjOz3PIQdMuHur4/1VD+vpXZWsVJyvKhru9PNZS/b2W2VvHlPjMzy62qTlKSDpL0qqRpks4qsr6dpFvT+mck9axAmGZmVatqL/dJag1cCewPzAQmSBodEVMKNvsu8GFEfFHSMcDPgaPLH601tY7rti06z5SZNV7Hdds2eZtVm6SA3YFpETEdQNItwOFAYZI6HLgwPb8D+J0kRUSUM1BrgMIBGMuXQqviFwuGbte1jEGZWWNVc5LqBrxd8HomMLCubSJimaT5QGfgg8KNJJ0InAiwxRZbNFe8VorCARgzJ8ICT25o1pJVc5JqMhFxNXA1wIABA9zLyouDf1bpCMxsDVXzwIlZQI+C193TsqLbSGoDdALmliU6MzOr6iQ1AdhGUi9J6wDHAKNrbTMaOD49Pwp4zPejzMzKp2ov96V7TD8EHgRaA3+KiMmSLgYmRsRo4DrgRknTgH+TJTIzMyuTqk1SABExBhhTa9n5Bc+XAF8vd1xmZpap5st9ZmaWc05SZmaWW05SZmaWW05SZmaWW/KI6qYlaQ7wZiN370KtahY54bgaJo9x5TEmcFwNtTbHtWVEbFJ7oZNUjkiaGBEDKh1HbY6rYfIYVx5jAsfVUNUYly/3mZlZbjlJmZlZbjlJ5cvVlQ6gDo6rYfIYVx5jAsfVUFUXl+9JmZlZbrknZWZmueUkZWZmueUklROSDpL0qqRpks6qUAw9JI2VNEXSZEmnpuUbS3pY0r/SvxtVKL7Wkp6XdG963UvSM+mc3ZqmXCl3TBtKukPSK5KmStozD+dL0unpZ/iypJslta/E+ZL0J0nvS3q5YFnR86PMFSm+FyXtUua4Lks/xxcl3S1pw4J1Z6e4XpV0YDnjKlj3X5JCUpf0uiznq66YJJ2cztdkSb8oWN605yoi/Kjwg2yqkNeBrYB1gBeA3hWIY3Ngl/S8I/Aa0Bv4BXBWWn4W8PMKnacfATcB96bXtwHHpOdXAd+vQEx/Bk5Iz9cBNqz0+QK6AW8A6xacpxGVOF/AYGAX4OWCZUXPD/Bl4H5AwB7AM2WO6wCgTXr+84K4eqffyXZAr/S72rpccaXlPcimFXoT6FLO81XHuRoKPAK0S6+7Nte5ck8qH3YHpkXE9Ij4FLgFOLzcQUTEuxHxXHq+EJhK9oF3ONmHMenfr5Y7NkndgUOAa9NrAV8C7qhUXJI6kf0CXwcQEZ9GxDxycL7IpuFZN80ovR7wLhU4XxHxONlcbIXqOj+HAzdE5mlgQ0mblyuuiHgoIpall0+TzdZdE9ctEfFJRLwBTCP7nS1LXMlvgP8GCke6leV81RHT94GfRcQnaZv3C2Jq0nPlJJUP3YC3C17PTMsqRlJPoD/wDLBpRLybVr0HbFqBkH5L9kv6WXrdGZhX8KFSiXPWC5gDXJ8uQ14raX0qfL4iYhbwS+AtsuQ0H3iWyp+vGnWdnzz9HnyHrJcCFY5L0uHArIh4odaqSsa1LbBPunz8d0m7NVdMTlK2CkkdgDuB0yJiQeG6yPr0Zf3egqRDgfcj4tlyHrcEbcgug/whIvoDH5FdvlqhQudrI7K/aHsBXwDWBw4qZwylqsT5qY+kc4FlwKgcxLIecA5wfn3bllkbYGOyy4xnArelqxtNzkkqH2aRXXOu0T0tKztJbckS1KiIuCstnl1zGSH9+35d+zeTQcBhkmaQXQr9EnA52eWNmtmlK3HOZgIzI+KZ9PoOsqRV6fO1H/BGRMyJiKXAXWTnsNLnq0Zd56fivweSRgCHAselBFrpuLYm+2PjhfT/vzvwnKTNKhzXTOCudKnxn2RXOLo0R0xOUvkwAdgmjb5aBzgGGF3uINJfQtcBUyPi1wWrRgPHp+fHA38tZ1wRcXZEdI+InmTn5rGIOA4YCxxVwbjeA96WtF1aNAyYQoXPF9llvj0krZd+pjVxVfR8Fajr/IwGvpVGre0BzC+4LNjsJB1Edkn5sIj4uFa8x0hqJ6kXsA3wz3LEFBEvRUTXiOiZ/v/PJBvc9B6VPV/3kA2eQNK2ZIOGPqA5zlVzjAbxo1EjaL5MNprudeDcCsWwN9mllxeBSenxZbL7P48C/yIb0bNxBc/TED4f3bdV+gWYBtxOGmlU5nj6ARPTObsH2CgP5wu4CHgFeBm4kWy0VdnPF3Az2X2xpWQfsN+t6/yQjVK7Mv0OvAQMKHNc08jup9T837+qYPtzU1yvAgeXM65a62fw+ei+spyvOs7VOsBf0v+v54AvNde5clkkMzPLLV/uMzOz3HKSMjOz3HKSMjOz3HKSMjOz3HKSMjOz3HKSMmsCyqqh/2fB6y9IumN1+zSg7QslnZGeXyxpvyZqd3OlivLNQdKiBmz7iCpUXd/yzUnKrGlsCKxIUhHxTkQcVffmjRMR50fEI03U3I+Aa5qorTV1IwXnz6yGk5RZ0/gZsLWkSWleop418+9IGiHpnjR30gxJP5T0o1SU9mlJG6fttpb0gKRnJT0hafvaB5E0UtJR6fkMSRdJek7SSzXbS1o/zQH0z3SMuirqHwk8kPa5T9LO6fnzks5Pzy+W9L30/ExJE9LcRRcVxPTNdKxJkv4oqXWtmLtIekrSIan39nja9mVJ+6TNRgPDG3nubS3mJGXWNM4CXo+IfhFxZpH1OwFHALsBlwAfR1aU9ingW2mbq4GTI2JX4Azg9yUc94OI2AX4Q9oHsm/8PxYRu5OVrrksVWdfIZWs+TDSVAvAE2RVrTuRFVcdlJbvAzwu6QCyEje7k1XZ2FXSYEk7AEcDgyKiH7AcOK7gOJsC9wHnR8R9wLHAg2nbvmSVHYiID4F2kjqX8J6tirSpfxMzawJjI5uja6Gk+cDf0vKXgJ1T5fm9gNsLikm3K6HdmiLAz5IlQcgm7zus5j4W0B7Ygmx+sBqbk00zUuMJ4BSyyRLvA/ZPFbh7RcSrqTd1APB82r4DWdLaGdgVmJDiXpfPC8a2JSt/9IOI+HtaNgH4UypkfE9ETCqI4X2yqu1zS3jfViWcpMzK45OC558VvP6M7PewFdl8T/0a2e5yPv99FnBkRLy6mv0WkyWvGhOAAcB04GGyitbfI0t+NW1eGhF/LGxE0snAnyPi7CLHWJb2PxD4O2QT6EkaTDaB5UhJv46IG9L27VNcZiv4cp9Z01gIdGzszpHN2/WGpK9DVpFeUt9GNvcgcHKqgI6k/kW2eQ3oWXD8T8mKq36d7BLkE2SXDx8vaPM7qceHpG6SupL1lI5Kz5G0saQta5olmzxwe0n/k9ZvCcyOiGvIZlnepeb9ApuRFVA1W8FJyqwJRMRcYHwaDHBZI5s5DviupBeAyWQTFzbGT8gutb0oaXJ6XTvej4DXJX2xYPETZJNLLk7Pu6d/iYiHgJuApyS9RDZ3VseImAKcBzwk6UWyXtjmBcdZTjYg4ktpiP4QsrmRnie7l3V52nRX4On4fOZgMwBXQTerVpK+BuwaEeflIJbLgdER8WilY7F88T0psyoVEXfnaDTdy05QVox7UmZmllu+J2VmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrnlJGVmZrlVsSQl6U+S3pf0cgP2GSIpJJ1QsKxfWnZGPfueJOlb9WzTT9KXS4hjgKQrSo27nrZGSPpdU7RlZra2qWRPaiRwUCP2exn4RsHr4cAL9e0UEVdFxA31bNYPqDdJRcTEiDilvu3MzGzNVCxJRcTjwL8bseubQHtJm0oSWaK7v2alpO9JmiDpBUl3SlovLb+wprclaZykn0v6p6TXJO0jaR3gYuBoSZMkHS1pd0lPSXpe0pOStkv7D5F0b0G7f0ptTpd0SkEs30zHmCTpj5Jap+XfTsf9JzCoMefPzKwa5OqelKQz0wd67UftS2t3AF8H9gKeAz4pWHdXROwWEX2BqcB36zhcm4jYHTgNuCAiPgXOB26NiH4RcSvwCrBPRPRP6/63jra2Bw4EdgcukNRW0g7A0cCgiOgHLAeOk7Q5cBFZctob6F3i6TEzqzptKh1AoYi4DLishE1vA24lSw43kyWrGjtJ+imwIdABeLCONu5K/z4L9Kxjm07AnyVtAwTQto7t7ouIT4BPJL0PbAoMA3YFJmQdPtYF3gcGAuMiYg6ApFuBbVfzXs3MqlaukpSkM4Hjiqx6vPAeUES8J2kpsD9wKisnqZHAVyPiBUkjgCF1HK6m97Wcus/DT4CxEfE1ST2BcfW0VdiegD9HxNmFG0r6ah1tmJlZLblKUg3oSUF2+a1rRCxPPZUaHYF3JbUlS3izGhDCwrR/jU4F+49oQDsAjwJ/lfSbiHhf0sap7WeAyyV1BhaQXbasd+CHmVk1quQQ9JuBp4DtJM2UVNe9o6Ii4smIuKfIqh+TJYLxZPeUGmIs0Ltm4ATwC+BSSc/TwIQeEVOA84CHJL0IPAxsHhHvAheSvffxZPfNzMysCEVEpWMwMzMrKlej+8zMzAo5SZmZWW616CQlqY2kOZJ+Vmv5OQ1o41pJdX5XKX1Jd0Aj4xss6TlJyyQdVWvdA5Lm1XwpOA8kbVfr+2kLJJ2Wg7japy9FvyBpsqSLKh0TNK60Vzk4roZxXKWrREwtOkmRDUF/Dfi6Vh7iV1KSktQ6Ik5Igxyaw1tkowJvKrLuMuA/mum4jRIRr6YvMvcj+47Xx8DdlY0KyIb4fyl9QbsfcJCkPSobEtD40l7NbSSOqyFG4rhKNZIyx9TSk9Rw4HKyZLAnQOpVrZt6AqNq7yBpkaRfSXoB2LOmpySptaSRkl6W9JKk02vt1yqt/2mpwUXEjIh4EfisyLpHyYa859Uw4PWIeLPSgURmUXrZNj0qPuJnDUp7NSvH1TCOq3SViClX35NqCEntgf2A/0dWXWI48GREnCXph6k3UMz6wDMR8V+pnZrl/YBuEbFTWr5hwT5tgFHAyxFxSZO+kfw6hqyaRy6kuofPAl8EroyIZyockpmVQUvuSR1KVg1iMXAn8NWaAq71WJ62r206sJWk/5N0ENkXbWv8kSpKUMqK7R4G3F7pWGpExPL0h0d3YHdJO1U4JDMrg5acpIYD+0maQfYXdmfgSyXstyQiltdeGBEfAn3JSh+dBFxbsPpJYGjqvVWDg4HnImJ2pQOpLSLmkX3pOm/X6s2sGbTIJCVpA2AfYIuI6BkRPYEfkCUugKWpLFJD2uwCtIqIO8kqRexSsPo6YAxwm6QWe4m0AYaTr0t9m9RcfpW0LtmAmYZWEzGzFqhFJinga8BjqfJ4jb8CX5HUDrgaeLHYwInV6AaMkzQJ+AuwUmHYiPg18Dxwo6SSzpuk3STNJKvP90dJkwvWPUF2OW1YKgt1YANibTaS1idLAnfVt20ZbQ6MTeWlJgAPR0TFh+6vaWmv5uK4GsZx5Tsml0UyM7Pcaqk9KTMzqwJOUmZmlltOUmZmlltlT1KSNpN0i6TXJT0raYyk3E6fLmlGGvnX1O2eLWmapFfzMmgCclsvrIeksZKmpNp9p1Y6Jsh1TUHH1QB5jCuPMUGF4oqIsj3IplR/CjipYFlfYJ8yx9GmAdvOALo08fF7k83G2w7oBbwOtC7nOVhNbIPJht+/XOlYCmLaHNglPe9IVq+xdw7iEtAhPW9LNtnmHo7Lca2NMVUqrnL3pIYCSyPiqpoFEfFCRDwBIOlMSRMkvViToSX1lDRV0jUpcz+UviuDpFPSX9cvSrolLdtY0j1p2dOSdk7LL5R0o6TxZMPIN5F0ZzreBEmD0nad0zEmS7qW7IfS1A4HbomITyLiDWAasHszHKfBIp/1wt6NiOfS84Vksxl3q2xUua4p6LgaII9x5TEmqExc5U5SO5FVh1iFpAOAbcg+rPsBu0oanFZvQ1avbUdgHnBkWn4W0D8idiarEgFwEfB8WnYOcEPBYXoD+0VETWHa30TEbqm9mgoTFwD/SMe6G9hiTd5wHboBbxe8nkkOPnRbAkk9gf5kf8FVnLLCxJOA98m+v+W4VsNxlS6PMUH548rTwIkD0uN54Dlge7LkBPBGRExKz58FeqbnLwKjJH0TWJaW7Q3cCBARjwGdlVWoABgdWa0/yIrT/i6d7NHABpI6kF3u+kva/z7gwyZ9l9Zo6edzJ3BaRCyob/tyiJzWFHRcDZPHuPIYE5Q/rnInqclk8xQVI+DSSPMZRcQXI+K6tK6wssRyPq/efghwJdk9lAmqv2TRRwXPW5FdS605XreCbmxzmwX0KHjdPS2zOigrc3UnMCoi8lQNA8hvTUHH1TB5jCuPMUH54ip3knoMaCfpxJoFknaWtA/wIPCd9NcykrpJ6lpXQ8pKE/WIiLHA/wCdgA7AE8BxaZshwAd1/NX9EHByQXv90tPHgWPTsoOBjRrzRusxGjhGUjtJvch6jP9shuOsFSSJrH7i1MjKU+WCclpT0HE1TB7jymNMKZayx1XWYqkREZK+BvxW0v8AS8hGz50WEf+StAPwVPaZxCLgm2Q9p2JaA3+R1ImsF3ZFRMyTdCHwJ2V13j4Gjq9j/1OAK9N2bciS00lk97RuVlZn70myCRWbVERMlnQbMIXsMuUPokhl9kpQVptrCNBFWd3BCwp6tJUyiGwW45fS5VmAcyJiTOVCArJRh39WNkVMK+C2yEFNQRxXQ+UxrjzGBBWIy7X7zMwst/I0cMLMzGwlTlJmZpZbTlJmZpZbLSJJSRqnrMbdpPS4o4nb7ynp2KZss4Rj5q52n3JaLwxA0oaS7pD0irIKJHvmIKbtCv5PTpK0QNJpjqvlxJXHmBxXrWO2hIETksYBZ0TExGZqf0hq/9DmaL/I8XqTTc++O/AF4BFg20qP8EtDvdePiEXpe0n/AE6NiKcrGReApD8DT0TEtZLWAdZL39PIhTTaaRYwMCLerHQ8NRxX6fIYEziuFtGTKkZSJ0lvpu9LIWl9SW9Laitpa0kPKKuy/oSk7dM2IyVdIelJSdMlHZWa+xmwT/rL4HRJO6YexSRlNQC3qSuORspl7b681gtLXzMYTPZdKSLi0zwlqGQY8HqePkQSx1W6PMYEVR5XS0pSowq6mJdFxHxgErBvWn8o8GBELAWuBk6OiF2BM4DfF7SzOVnppEPJkhNkNQCfSJUnfkP2fanLU+mPAWS19ZpSbmv3KZ/1wnoBc4DrJT0v6VpJ61c6qFqOIesd543jKl0eY4Iqj6slJanjCkoYnZmW3QocnZ4fA9yqrGLFXsDt6cP2j2SJqcY9EfFZREwBNq3jWE8B5yj7wvGWBfX+1no5rRfWhqz01R8ioj9ZeauzKhvS59Llx8OA2ysdSyHHVbo8xgSOC1pWkipmNHCQpI3JagI+Rvae5hUktH4RsUPBPoV1AItOwxERN5H9ABYDYyR9qYnjzn3tvpzVC5sJzCzo1d1BlrTy4mDguYiYXelAanFcpctjTOC4WnaSSvdPJpBNu3Fv6gUsAN6Q9HXIBgNI6ltPUwvJJtMj7bMVMD0irgD+CuzcxKHnsnafclovLCLeA96WtF1aNIyspFReDCefl2McV+nyGBM4rhY1um9zsp4NZEVj90vrjiLrcg6JiL+nZb2AP6R92pINUrhY0kiyZHZH2m5RRHRII9keBDoDI8lmzP0PYCnwHnBsRDTpRICSzgW+Q1a777SIuL8p228MZRNE/pmsLmJNXa6LKxtVRlkB4GuBdYDpwLcjouLTqKR7Y28BW6X7pLnguEqXx5jAca04XktIUmZmVp1a9OU+MzNbuzlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbjlJmZlZbq02SUnqIWmspCmSJks6tdSGJQ2RFJK+UrDsXklD6tnvNEnrlXqchpD0ZAnblHR8SddK6t1Ecc2Q1KUp2jIzW5vU15NaBvxXRPQG9gB+0MAP5pnAuQ2M6TSgWZJUROzVVMePiBMiIk+zw5qZrXVWm6Qi4t2IeC49XwhMBbo1oP0XgPmS9q+9QtIwSc9LeknSn9JU6qcAXwDGShpbZJ8Zki6VNEnSREm7SHpQ0uuSTkrbdJD0qKTnUtuHF+y/KP07RNI4SXdIekXSqDTN/CrHl/SHdKzJki4qaGucpAE17Uq6RNILkp6WtGlavomkOyVNSI9BaXlnSQ+lNq8F1IBzamZWPSKipAfQk2zK4A3S6zOBSUUeV6T1Q4B7gcHA39Oye9Py9sDbwLZp+Q1kU6gDzAC61BHDDOD76flvgBeBjsAmwOy0vE1BjF2AaXw+A/GigtjmA93JEvVTwN7Fjg9snP5tDYwDdk6vxwED0vMAvpKe/wI4Lz2/qaDdLYCp6fkVwPnp+SFp/6Lv2Q8//PCjmh9t6slhQNY7Ae5MiWQBQERcBlxW374R8bgkJO1dsHg74I2IeC29/jPwA+C3JYQzOv37EtAhsh7eQkmfSNoQ+Aj4X0mDgc/Ien6bAu/VauefETEzvb9JZEn4H0WO9w1JJ5Ilv82B3mTJsdCnZAkY4Fmgpue4H9BbWtFR2iCdy8HAEQARcZ+kD0t432ZmVafeJCWpLVmCGhURdxUsPxM4rsguj0fEKbWWXQKcR3aPa019kv79rOB5zes2KaZNgF0jYqmkGWQ9t7raAVhOkXMhqRdwBrBbRHwoaWQdbS2NiCjSVitgj4hYUqvdOt+cmZl9rr7RfQKuI7tM9evCdRFxWUT0K/KonaCIiIeAjYCd06JXgZ6Svphe/wfw9/R8IdklvMbqBLyfEtRQYMsG7l94/A3Iembz032mgxvY1kPAyTUvJPVLTx8Hjk3LDiY7N2ZmVkt9o/sGkSWQL6XBCpMkfbmRx7oE6AGQehbfBm6X9BJZL+iqtN3VwAPFBk6UaBQwILX7LeCVBu6/4vgR8QLwfGrjJmB8A9s6JcXyoqQpwElp+UXAYEmTyS77vdXAds3MqoI+v0plZmaWL644YWZmueUkZWZmudWik5SkNpLmSPpZreXnNKCN1ZY3KvzSbiPiG5y+VLxM0lEFy/tJeip9mfdFSUc3pv3mIGnDgi85T5W0Zw5ianR5ruYm6SBJr0qaJumsSsdTw3GVLo8xgeNaodJf1FqTB9lou/HA66T7a2n5ohL3b13CNuNIX9ptRHw9yUY03gAcVbB8W2Cb9PwLwLvAhpU+nymePwMnpOfr5CEusu+n7ZKedwReA3rnIK7W6f/eVulcveC4WlZceYzJca38aNE9KWA4cDnZ6Lg9AVKvat00EnFU7R1SCaNfSXoB2LOmpySptaSRkl5WVk7p9Fr7tUrrf1pqcBExIyJeJBu9WLj8tYj4V3r+DvA+2Xe7KkpSJ7IvGl8HEBGfRsS8igZFk5Tnai67A9MiYnpEfArcAhxezz7l4LhadkzguFZosUlKUnuyig5/A24mS1hExFnA4si+s1Xsy8brA89ERN+IKKww0Q/oFhE7RUQf4PqCdW3Ihrb/KyLOa+L3sTvZXySvN2W7jdQLmANcr6yu4rWS1q90UIUk9QT6A89UOBTIEuXbBa9nko/k6bhKl8eYwHGt0GKTFHAoMDYiFpNVxPiqpNYl7Lc8bV/bdGArSf8n6SBgQcG6PwIvR8Qlaxp0IUmbAzcC346Iz+rbvgzaALsAf4iI/mRfZM7TtfBVynOZ2dqtJSep4cB+qezRs0Bn4Esl7LckIpbXXhgRHwJ9ye5BnQRcW7D6SWBo6r01CUkbAPcB50bE003V7hqaCcyMiJpeyh1kSavi6irPVWGzSF9QT7qnZZXmuEqXx5jAca3QIpNU+oDfB9giInpGRE+yArXD0yZL04daQ9rsArSKiDvJ6gwWfjhfB4wBbpNUUlHeeo61DnA3cENE3LGm7TWViHgPeFvSdmnRMKDic2atrjxXhU0AtpHUK/1Mj+HzAsiV5LhadkzguFZY4w/cCvka8FhEFBaJ/SvwC0ntyEobvSjpuTruSxXTjexeTE3iPrtwZUT8Og0suFHScaVcnpO0G1ky2gj4iqSLImJH4BtkAxQ6SxqRNh8REZNKjLU5nQyMSv8Bp5OVr6q0mvJcLymrWA9wTkSMqVxIEBHLJP0QeJBs1NOfImJyJWMCx9XSYwLHVchlkczMLLda5OU+MzOrDk5SZmaWW05SZmaWW2VPUpI2k3SLpNclPStpjKRtyx1HqSTNSCP/mrrds1Ptq1clHdjU7TeWpD9Jel/Sy5WOpVAe48pjTOC4Gspxla4SMZU1SaWhxHcD4yJi64jYlWwU3aZljqOioxqVFbQ9BtgROAj4fYlfRC6HkWQx5c1I8hfXSPIXEziuhhqJ4yrVSMocU7l7UkOBpRFRMwsvEfFCRDwBIOlMSROUVQa/KC3rqawa9zXKKmA/JGndtO4UZZWxX5R0S1q2saR70rKnJe2cll8o6UZJ48mGkW8i6c50vAmSBqXtOqdjTJZ0LaBmOA+HA7dExCcR8QYwjawmVsVFxOPAvysdR215jCuPMYHjaijHVbpKxFTuJLUTWXWIVUg6ANiG7MO6H7CrpMFp9TbAlek7RvOAI9Pys4D+EbEzK0/N/nxadg5ZBfIavYH9IqKmMO1vImK31F5NhYkLgH+kY90NbLEmb7gOea3LZWaWK3n6Mu8B6fF8et2BLDm9BbxR8EXXZ8mmwAB4keyLp/cA96Rle5OSWEQ8lnpGG6R1o1OtP8iK0/bOrkACsEGqDTcYOCLtf5+kD5vuLZqZWUOUO0lNBo6qY52ASyPijystzKpeF1aWWA6sm54fQpZUvgKcK6lPPcf/qOB5K2CPiFhS63j1NNEk8lqXy8wsV8p9ue8xoJ2kE2sWSNpZ0j5kZTa+k3ozSOomqWtdDaXyRT0iYizwP0Anst7XE8BxaZshwAd1VMx+iKwEUE17/dLTx4Fj07KDyUoaNbXRwDGS2knqRdZj/GczHMfMrEUra5KKrAbT18iql78uaTJwKfBeRDwE3AQ8JeklsgrcHVfTXGvgL2nb54Er0gR9F5Ldz3oR+BlwfB37nwIMSAMsprDyPa3BKbYjyC43NqlU6+o2suKtDwA/KFaZvRIk3Qw8BWwnaaak71Y6JshnXHmMCRxXQzmufMfk2n1mZpZbrjhhZma55SRlZma55SRlZma51SKSlKRxqcbdpPRo0tlsU1WLY5uyzRKOmdfafTMkvZTO88RKx1ND0qmSXk6VQE4r87FXqVcm6espls8kDShnPI5r7YjJcZWmRSSp5LiI6JcedX3XqrF6koadl4PyXbsPYGg6zxX5BalN0k7A98iqkfQFDpX0xTKGMJJV65W9TDb68/EyxlHbSBxXqUaSv5jAcdWrJSWplUjqJOnN9H0pJK0v6W1JbSVtLekBZVXWn5C0fdpmpKQrJD0pabqkmmT3M2Cf1Hs4XdKOkv6ZXr8oaZsmDj+3tftyagfgmYj4OCKWAX8nVQUph2L1yiJiakS8Wq4YinFcpctjTCkGx1WPlpSkRhVc7rssIuYDk4B90/pDgQcjYilwNXByqrJ+BvD7gnY2JyuddChZcoKsBuATqffwG7LvTF0eEf2AAWS19ZpSnmv3BfBQSvAn1rt1ebxM9kdEZ0nrAV9m5YodZraWylPtvvocFxG175HcChwNjCW7fPZ7ZRUr9gJu1+cljtoV7HNPRHwGTJFU1xQhT5GVWeoO3BUR/2qqN9EC7B0Rs5RV+3hY0ivpr6qKiYipkn5OViXkI7I/TnLx5Wcza14tqSdVzGjgIEkbA7uSlV1qBcwruH/VLyJ2KNinsA5g0UJ9EXETcBiwGBgj6UtNHHdua/dFxKz07/tkVeBzcRkyIq6LiF0jYjDwIfBapWMys+bXopNURCwCJpBNu3FvRCxPdfrekPR1yCZalNS3nqYWUlCCSdJWwPSIuAL4K7BzE4eey9p96b5ex5rnZFXpczEraOrZIWkLsvtRN1U2IjMri4jI/QMYB7xKdplnEvBIwbqjyO6j7FuwrBdZTbwXyOrjnZ+WjwSOKthuUfq3LVkv7AXgdLJ7VJPTsR4ANm6G93Qu8Hp6XwdX+hynmLZK5+CF9P7PrXRMBbE9kX6WLwDDynzsm4F3gaVk9w+/S1aDciZZz3w22f3Qcp8Tx9WCY3JcpT1cu8/MzHKrRV/uMzOztZuTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5ZaTlJmZ5db/B2Lq2SvCd5uvAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "T1, T2 = df['time']  [df['group'] == 1], df['time']  [df['group'] == 2]\n",
    "E1, E2 = df['status'][df['group'] == 1], df['status'][df['group'] == 2]\n",
    "kmf1, kmf2 = fit_plot(T1, T2, E1, E2, 'Kaplan-Meier survival curves', 'weeks', '1=Maintained', '2=Not maintained')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The vertical ticks on the graphs indicate censored events. At the bottom, the number of patients \"at risk\" is shown. The study concerned 11 patients in group 1 with one patient censored after 161 weeks and 12 patients in group 2. A Kaplan-Meier curve estimates the survival probability as a function of time (duration). \n",
    "\n",
    "In this example, the two curves do not appear to differ very much. To analyze this more precisely one performs a [logrank test](https://en.wikipedia.org/wiki/Logrank_test):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.06533932204050484"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lifelines.statistics.logrank_test(T1, T2, E1, E2).p_value"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The null hypothesis is that survival is the same for both groups. Since $p=0.065$ is not particularly small (e.g., not below $\\alpha=0.05$), the null hypothesis is not strongly rejected. In other words, the logrank test also says that the curves do not differ significantly---with the usual caveat of small sample sizes."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Privacy-Preserving Survival Analysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In a multiparty setting, each party holds a private dataset and the goal is to perform a survival analysis on the *union* of all the datasets. To obtain the secure union of the datasets in an efficient way, each private dataset will be represented using two survival tables, one survival table per group of patients. \n",
    "\n",
    "The survival tables for the aml dataset are available from the Kaplan-Meier fitters `kmf1` and `kmf2` output by the call to `fit_plot()` above:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>removed</th>\n",
       "      <th>observed</th>\n",
       "      <th>censored</th>\n",
       "      <th>entrance</th>\n",
       "      <th>at_risk</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>event_at</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0.0</th>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>11</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>13.0</th>\n",
       "      <td>2</td>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>18.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>28.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>31.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>34.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>45.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>48.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>161.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          removed  observed  censored  entrance  at_risk\n",
       "event_at                                                \n",
       "0.0             0         0         0        11       11\n",
       "9.0             1         1         0         0       11\n",
       "13.0            2         1         1         0       10\n",
       "18.0            1         1         0         0        8\n",
       "23.0            1         1         0         0        7\n",
       "28.0            1         0         1         0        6\n",
       "31.0            1         1         0         0        5\n",
       "34.0            1         1         0         0        4\n",
       "45.0            1         0         1         0        3\n",
       "48.0            1         1         0         0        2\n",
       "161.0           1         0         1         0        1"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>removed</th>\n",
       "      <th>observed</th>\n",
       "      <th>censored</th>\n",
       "      <th>entrance</th>\n",
       "      <th>at_risk</th>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>event_at</th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "      <th></th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0.0</th>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5.0</th>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8.0</th>\n",
       "      <td>2</td>\n",
       "      <td>2</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>12.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>16.0</th>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>7</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>23.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>6</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>27.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>5</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>30.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>4</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>33.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>3</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>43.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>2</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>45.0</th>\n",
       "      <td>1</td>\n",
       "      <td>1</td>\n",
       "      <td>0</td>\n",
       "      <td>0</td>\n",
       "      <td>1</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "          removed  observed  censored  entrance  at_risk\n",
       "event_at                                                \n",
       "0.0             0         0         0        12       12\n",
       "5.0             2         2         0         0       12\n",
       "8.0             2         2         0         0       10\n",
       "12.0            1         1         0         0        8\n",
       "16.0            1         0         1         0        7\n",
       "23.0            1         1         0         0        6\n",
       "27.0            1         1         0         0        5\n",
       "30.0            1         1         0         0        4\n",
       "33.0            1         1         0         0        3\n",
       "43.0            1         1         0         0        2\n",
       "45.0            1         1         0         0        1"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "display(kmf1.event_table, kmf2.event_table)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Using `mpc.input()`, each party will secret-share its pair of survival tables with all the other parties. To allow for a simple union of the survival tables, however, we modify the representation of the survival tables as follows. \n",
    "\n",
    "First, we determine when the last event happens, which is at $t=161$ in the example. We compute $maxT$ securely as the maximum over all parties, using the following MPyC one-liner:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "161"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "secfxp = mpc.SecFxp(64)  # logrank tests will require fixed-point arithmetic\n",
    "maxT = int(await mpc.output(mpc.max(mpc.input(secfxp(int(df['time'].max()))))))\n",
    "timeline = range(1, maxT+1)\n",
    "max(timeline)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "All parties use $1..maxT$ as the global timeline. Subsequently, each party pads its own survival tables to cover the entire timeline. This is accomplished by two calls to the function `events_to_table()`, which only keeps the essential information:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>d1</th>\n",
       "      <th>n1</th>\n",
       "      <th>d2</th>\n",
       "      <th>n2</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>5</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>2</td>\n",
       "      <td>12</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>6</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>7</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>8</th>\n",
       "      <td>0</td>\n",
       "      <td>11</td>\n",
       "      <td>2</td>\n",
       "      <td>10</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>9</th>\n",
       "      <td>1</td>\n",
       "      <td>11</td>\n",
       "      <td>0</td>\n",
       "      <td>8</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>10</th>\n",
       "      <td>0</td>\n",
       "      <td>10</td>\n",
       "      <td>0</td>\n",
       "      <td>8</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "    d1  n1  d2  n2\n",
       "1    0  11   0  12\n",
       "2    0  11   0  12\n",
       "3    0  11   0  12\n",
       "4    0  11   0  12\n",
       "5    0  11   2  12\n",
       "6    0  11   0  10\n",
       "7    0  11   0  10\n",
       "8    0  11   2  10\n",
       "9    1  11   0   8\n",
       "10   0  10   0   8"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d1, n1 = events_to_table(maxT, T1, E1)\n",
    "d2, n2 = events_to_table(maxT, T2, E2)\n",
    "pd.DataFrame({'d1': d1, 'n1': n1, 'd2': d2, 'n2': n2}, index=timeline).head(10)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Column `d1` records the number of observed events (\"deaths\") for group 1 on the entire timeline, and column `n1` records the number of patients \"at risk\". Similarly, for group 2. \n",
    "\n",
    "To obtain the secure union of the privately held datasets, we let all parties secret-share their survival tables and simply add all of these elementwise:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "d1, n1, d2, n2 = (functools.reduce(mpc.vector_add, mpc.input(list(map(secfxp, _)))) for _ in (d1, n1, d2, n2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The joint dataset (union) is now represented by `d1, n1, d2, n2`. Note that these values are all secret-shared, for example:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<mpyc.sectypes.SecFxp64:32 at 0x26f1269f820>"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "d1[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Aggregate Kaplan-Meier Curves\n",
    "\n",
    "We now proceed to analyze the joint dataset in a privacy-preserving way by plotting *aggregated* versions of the Kaplan-Meier curves. The exact Kaplan-Meier curves would reveal too much information about the patients in the study. By aggregating the events over longer time intervals, the amount of information revealed by the curves is reduced. At the same time, however, the aggregated curves may still be helpful to see the overall results for the study---and in any case to serve as a sanity check.\n",
    "\n",
    "The function `aggregate()` securely adds all the observed (\"death\") events over intervals of a given length (stride). The aggregated values are all output publicly, and used to plot the curves via the function `events_from_table()`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAakAAAGcCAYAAABjm19FAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAABDfklEQVR4nO3deZwU1bn/8c+XRXBBjAKKoIIRF1wARXEPxD0xmkSNC7nKTYw/s7jlauIWt2g012x6Y2JcItHgmqghirsQDEoCKC6ACyoKiopGBVSU5fn9UWewGHqYGZzprqa/79erX1N1quqcp6tn+plTyylFBGZmZkXUptIBmJmZNcRJyszMCstJyszMCstJyszMCstJyszMCstJyszMCstJyiyRNEbSsRVod4akvcvdblNJGirp/krH0VySrpT0kxaoZ7ikC1siJms+JykDln5BvyupQ6VjWRmSekkKSe1aqf7zJP05N99D0rOSLpek1mizJaR98lZ+v0hqn8qadJNkRIyIiH1bL8rWERHHR8RPKx2HfTZOUoakXsAeQAAHtWI7rZJAyk3SJsBYYGREnBjFvyP+XeCA3PwBqazVSWrbSvVK0irz/bWq/G20hlXmQ7bP5GhgPDAcOCa/QNJ6kv4uaa6kCZIulPTP3PJ9JT0n6X1Jv5P0j7pDZpKGSRon6deS3gHOk9RB0i8kvSrpzXRIZvVcfT+SNFvS65KOTT2BzdKyL0t6IsUyU9J5uVDHpp/vSZovaZe0zbckTUu9xPtSgqlra5/UG3pf0m+BRntEkj6f2hoRET/KlV+WYporaZKkPXLLzpP0F0m3SJon6XFJ/RqofydJj0l6L+2H30paLbc8JB0v6YW0zhVN6MndQPYZ1zkauL5eu50lXZvafC19zm3TsmH1PvMtJT0g6T/ps/9GbtlwSb+XNErSB8CQEu9xmKSX0r54WdLQ3H7K91aX6R0r6+1fJGkc8CFwmqSJ9eo+RdLIXCwXpulpkg7MrddO0hxJ26f52yS9kX4XxkraupF9mm/zO6n+eZKm5upc+rtbIp7BkmZJ+rGkN4DrmhDjzpIeTZ/7k5IGN7ZPVwkR4VeNv4DpwPeAHYCFwPq5ZTen1xpAX2Am8M+0rAswF/g60A44KW1/bFo+DFgEnJCWrw78GhgJrAt0Av4OXJzW3x94A9g6tfdnst7dZmn5YGBbsn+utgPeBL6alvVK67bLxX5wem9bpfbPBh7NxT4POBRoD5ySYj22gX10HvAo8BpwZonl3wTWS+38T3ofHXPbLsy1dSrwMtA+LZ8B7J2mdwB2TvX0AqYBJ+faCeAuYB1gY2AOsP8KPtsAtkn7ah3gc2l6m+zPf+l6dwB/ANYEugH/Bv5f7nOs+8zXTL8D/51iHAC8DfRNy4cD7wO7pc+pY7141iT7ndkizXcHts7tpz/n1l3mMwXGAK+S/X60Azqnz7BPbpsJwBG5WC5M0+eQ/WNRt96XgWm5+W+R/T52AH4DTM4tW1pPif17WPqd2JHsn5zNgE1y+36zUvWQ/S4vAn6e2lx9RTECPYB3gC+l/bpPmu+6on26KrwqHoBfFf4FgN3JvkC7pPlngVPSdNu0bIvc+hfmvrCOBh7LLVP6AssnqVfrLf8A+HyubBfg5TT9R1LCSvOb1f9Drxf7b4Bfp+llvtBS2T3At3Pzbcj+A98kxT6+XmyzWHGSmgu8l49/Bfv1XaBfbtvx9eKYDeyR5meQklSJek4G7sjNB7B7bv5W4PQVxBFpP14D/D/geODqun2b1lkf+BhYPbfdkcDo3OdY95kfDjxSr40/AOem6eHA9SuIZ820Dw/Jt5fbT40lqQvqbfNn4Jw03Ycsaa2Ri6UuKWxWb9mIuu1KxLhOardz/XpKrHsfcNKK9n1uPh/PYOATckl8RTECPwZuKNH2MSvap6vCy4f77Bjg/oh4O83fyKeH/LqS/cc6M7d+fnrD/Hz6xptVr/78+l3JekiT0iGL94B7U/ly9dWbRtIgSaPTIZD3yb5wu6zgvW0CXJZr6z9kyahHA7HPLFVJzkiyRPpw/rBhiu3UdLjm/dRW53qx5dtaQrafNqzfgKTNJd2VDj3NBX5W4j2+kZv+EFgrbTtF2aHO+fnDjcn1ZIl5uUN9ZPupPTA7t6/+QNajqm8TYFDdemndocAGpd5rfRHxAVmiOz61d7ekLRtav4T6dd9IllABjgLujIgPS7Q7naxX+hVJa5Cde70RsvNmki6R9GLa5zPSZiv63aqzEfBiM+LPmxMRC5oSI9l+P6zeft8d6N4C+7TQfLKuhik7F/QNoG06Lg7ZoYd1lJ0zeYbskERP4Pm0fKNcFbPTsrr6lJ9P8hcVvA18RHYo4rUSIS1TX722IPuD/S1wQEQskPQbPv0iCZY3E7goIkbUXyCpT77+FHv99pYTET9UdgXkw5L2jIjXUkL4EbAXMCUilkh6l2XPceXbapPe5+slmvg98ARwZETMk3Qy2WHCRkXEis6jPEJ2GCiAfwKfzy2bSdaT6hIRixppZibwj4jYZ0WhNBLnfcB96ffvQrKe3R5kvew1cqtuUGrzevMPAF0l9SdLVqesoOmb0jptgKkpKUCW3A4G9iZLUJ3JesJNuWpzJsvuy7wPWf795P+JK7WfGopxJllP6julGlrBPq167knVtq8Ci8nONfVPr63IvtCOjojFwO1kFzyskf47y5+AvxvYVtJX08nt71P6iwVY2oO4Gvi1pG6w9FLu/dIqtwL/LWmr9J9k/XtcOgH/SQlqJ7IvlzpzgCXAprmyK4Ez6k6CK7s44LBc7FtL+nqK/cQVxV7PD4DRwEOS1k9xLUoxtJN0DrB2vW12yLV1MllSGF+i7k5khxXnp/393SbGtEKpp/gV4KA0nV82G7gf+KWktSW1kfR5SV8oUdVdwOaS/kvZpeztJe0oaaumxCFpfUkHS1qTbB/MJ/vcACYDe0raWFJn4IwmvK+FwG3ApWTnOR9Yweo3A/uS7dMbc+WdUizvkCWVnzXlvSTXAKdK2kGZzXK97MnAUamntj9Qan82NcY/k/Ww9kv1dUwXX/RsZJ9WPSep2nYMcF1EvBoRb9S9yHorQ9MX6g/I/rN8g+wqsZvI/hBIhwgPA/6X7A+8LzCxbnkDfkx2McP4dGjlQWCLVN89wOVkCWA6n36J19X3PeACSfPITjLfWldpOsRzETAuHQ7ZOSLuIDsxfXNq6xnSpdi52C9JsfcBxjVlp6Uv+ePILi54EJhEdtjyeeAVYAHLH5b6G9khmXeB/wK+nr5g6zuVLPnOI0votzQlpibGPSUipjSw+GhgNWBqivEvZD2v+nXMI/sSPYKsJ/gGn578b4o2wA/Ttv8h++L+bqr7AbL3+xTZPr2riXXeSNYLum1FPcGUjB8DdmXZ/Xo92ef2Gtn7L/XPQ0N13kb2e3cj2Wd2J1myhOxCoq+QnS8ampY1Vl/JGCNiJllv70yyf4ZmAqeR7c8G9+mqQPX+qTJbIUk/BzaIiGNKLGtDdjhjaESMboG2tiJLLB2acBiqsJRdKr9ZRHyz0rGYVRv3pGyFlN0Ts106lLET8G2yy5Xrlu8naZ10nuZMsuP4Tf5PtER7X1N2L9XnyP5D/3s1Jygz+2ycpKwxncjOS31Advjhl2SHrursQnZ109tkhza+GhEffYb2/h/wVqpzMavQYQszaz4f7jMzs8JyT8rMzArLScrMzArLN/O2sC5dukSvXr0qHYaZWVWZNGnS2xHRtX65k1QL69WrFxMnTmx8RTMzW0rSK6XKfbjPzMwKy0nKzMwKy0nKzMwKy+ekzKzsFi5cyKxZs1iwYEHjK9sqpWPHjvTs2ZP27ds3af2aTVKS/ggcCLwVEduUWC7gMrInYX4IDIuIx8sbpdmqadasWXTq1IlevXqR/alZLYgI3nnnHWbNmkXv3r2btE0tH+4bTva48oYcQDYydh+yEa9/X4aYzGrCggULWG+99Zygaowk1ltvvWb1oGs2SUXEWLJh7RtyMNljsCMixpM9CHC5RxeY2cpxgqpNzf3ca/ZwXxP0YNlnAs1KZbNbo7Hxv/sOnd6b1hpVN2h+n68x6LD/KWubZkXxrW99i7vuuotu3brxzDPPNGmbMWPGMGTIEK6++mqOPfZYACZPnsyAAQO49NJLOfXUUxvc9sorr2SNNdbg6KOPbnCdyZMn8/rrr/OlL31phXFMnDiR66+/nssvv7xJca/I8OHDmThxIr/97W8/c12toWZ7Ui1J0nGSJkqaOGfOnJWuZ8miRWV7bfTxi6wx7bYW3Atm1WXYsGHce++9zd5um2224dZblz5vk5tuuol+/fo1ut3xxx+/wgQFWZIaNWpUo3UNHDiwRRJUNXBPqmGvARvl5numsuVExFXAVQADBw5cqWHld/7e1Suz2Uqb8rPdYZEf02S1a88992TGjBnN3m6TTTZh7ty5vPnmm3Tr1o177713mZ7P1VdfzVVXXcUnn3zCZpttxg033MAaa6zBeeedx1prrcWpp57K4MGDGTRoEKNHj+a9997j2muvZdCgQZxzzjl89NFH/POf/+SMM86gd+/enHTSSSxYsIDVV1+d6667ji222IIxY8bwi1/8grvuuovzzjuPV199lZdeeolXX32Vk08+mRNPPBGAP//5z1x++eV88sknDBo0iN/97ne0bduW6667josvvph11lmHfv360aFDUx+sXH5OUg0bCfxA0s3AIOD99GhnM2tB5/99ClNfn9uidfbdcG3O/crWzd7u0ksvZcSIEcuV77nnnsv0XA499FBuu+02BgwYwPbbb7/Ml/zXv/51vvOd7wBw9tlnc+2113LCCScsV+eiRYv497//zahRozj//PN58MEHueCCC5Y59DZ37lweeeQR2rVrx4MPPsiZZ57JX//61+XqevbZZxk9ejTz5s1jiy224Lvf/S7Tp0/nlltuYdy4cbRv357vfe97jBgxgn322Ydzzz2XSZMm0blzZ4YMGcKAAQOava/KpWaTlKSbgMFAF0mzgHOB9gARcSUwiuzy8+lkl6D/d2UiNbNyOe200zjttNMaXe8b3/gGhx9+OM8++yxHHnkkjz766NJlzzzzDGeffTbvvfce8+fPZ7/99itZx9e//nUAdthhhwZ7dO+//z7HHHMML7zwApJYuHBhyfW+/OUv06FDBzp06EC3bt148803eeihh5g0aRI77rgjAB999BHdunXjX//6F4MHD6Zr12ws18MPP5znn3++0fdcKTWbpCLiyEaWB/D9MoVjVrNWpsfTWprak9pggw1o3749DzzwAJdddtkySWrYsGHceeed9OvXj+HDhzNmzJiSbdX1vtq2bcuiBg69/+QnP2HIkCHccccdzJgxg8GDB6+wrnx9EcExxxzDxRdfvMy6d955Z8k6iqpmk5SZWX1N7UkBXHDBBbz11lu0bdt2mfJ58+bRvXt3Fi5cyIgRI+jRo0eT2+/UqRPz5s1bOv/+++8v3X748OFNrgdgr7324uCDD+aUU06hW7du/Oc//2HevHkMGjSIk046iXfeeYe1116b2267rUkXflSKr+4zs5p05JFHsssuu/Dcc8/Rs2dPrr322mZtv+uuu/LVr351ufKf/vSnDBo0iN12240tt9yyWXUOGTKEqVOn0r9/f2655RZ+9KMfccYZZzBgwIAGe1sN6du3LxdeeCH77rsv2223Hfvssw+zZ8+me/funHfeeeyyyy7stttubLXVVs2qt9yUHdWyljJw4MCohudJTfnZ7ixZtIhtzxlf6VCsBk2bNq3wX47Wekp9/pImRcTA+uu6J2VmZoXlc1I1LICRk0ve+tUqOq3eniFbdCtbe2ZW/ZykalkEXTt1LFtzc+b5sQxm1jw+3GdmZoXlJGVmZoXlJGVmZoXlJGVmNWfmzJkMGTKEvn37svXWW3PZZZc1edsxY8Ygib///e9Lyw488MAGR5ao85vf/IYPP/xwZUNeoV133bXRdZra/rHHHsvUqVNbIix69erF22+//ZnqcJIys5rTrl07fvnLXzJ16lTGjx/PFVdc0awv5p49e3LRRRc1q83WTFL5YZk+a/vXXHMNffv2bYmwWoSTlJnVnO7du7P99tsD2VBEW221Fa+91vTbMfr160fnzp154IEHllv20EMPMWDAALbddlu+9a1v8fHHH3P55Zfz+uuvM2TIEIYMGbLcNr169eKMM86gf//+DBw4kMcff5z99tuPz3/+81x55ZUAzJ8/n7322ovtt9+ebbfdlr/97W9Lt19rrbWArJc3ePBgDj30ULbcckuGDh1KRJRs/7vf/S4DBw5k66235txzz11a1+DBg6kbkGCttdbirLPOol+/fuy88868+eabAMyZM4dDDjmEHXfckR133JFx48YB8M4777Dvvvuy9dZbc+yxx9ISg0X4EnQzq6x7Toc3nm7ZOjfYFg64pEmrzpgxgyeeeIJBgwYBTR9k9qyzzuInP/kJ++yzz9KyBQsWMGzYMB566CE233xzjj76aH7/+99z8skn86tf/YrRo0fTpUuXknFsvPHGTJ48mVNOOYVhw4Yxbtw4FixYwDbbbMPxxx9Px44dueOOO1h77bV5++232XnnnTnooIOWexz7E088wZQpU9hwww3ZbbfdGDduHCeeeOJy7V900UWsu+66LF68mL322ounnnqK7bbbbpm6PvjgA3beeWcuuugifvSjH3H11Vdz9tlnc9JJJ3HKKaew++678+qrr7Lffvsxbdo0zj//fHbffXfOOecc7r777mYPNVWKk5SZ1az58+dzyCGH8Jvf/Ia1114baPogs3vuuScA//znP5eWPffcc/Tu3ZvNN98cgGOOOYYrrriCk08+udH6DjroIAC23XZb5s+fT6dOnejUqRMdOnTgvffeY8011+TMM89k7NixtGnThtdee40333yTDTbYYJl6dtppJ3r27AlA//79mTFjBrvvvvty7d16661cddVVLFq0iNmzZzN16tTlktRqq63GgQceCGSPFKnrOT744IPLHB6dO3cu8+fPZ+zYsdx+++1A9viQz33uc42+78Y4SZlZZTWxx9PSFi5cyCGHHMLQoUOXPtsJmt6Tgqw3deGFF9Ku3Wf/Kq173EabNm2WefRGmzZtWLRoESNGjGDOnDlMmjSJ9u3b06tXLxYsWP4G+VKP7ajv5Zdf5he/+AUTJkzgc5/7HMOGDStZV/v27Zf21PJ1LVmyhPHjx9OxY+sPBuBzUmZWcyKCb3/722y11Vb88Ic/XGbZaaedxuTJk5d71U9QAPvuuy/vvvsuTz31FABbbLEFM2bMYPr06QDccMMNfOELXwCWfwxHc73//vt069aN9u3bM3r0aF555ZVmbZ9vf+7cuay55pp07tyZN998k3vuuadZde2777783//939L5yZMnA1kiv/HGGwG45557ePfdd5tVbylOUmZWc8aNG8cNN9zAww8/TP/+/enfvz+jRo1aqbrOOussZs6cCUDHjh257rrrOOyww9h2221p06YNxx9/PADHHXcc+++/f8kLJ5pi6NChTJw4kW233Zbrr7++2Y8Bybffr18/BgwYwJZbbslRRx3Fbrvt1qy6Lr/8ciZOnMh2221H3759l17cce655zJ27Fi23nprbr/9djbeeONm1VuKH9XRwqrpUR2LFy7kg2+u3B/mypgzbwEH9W/6A+Bs1eVHddQ2P6rDzMxWCU5SZmZWWE5SZmZWWE5SZlYRPh9em5r7uTtJmVnZdezYkXfeeceJqsZEBO+8806z7q/yzbxmVnY9e/Zk1qxZzJkzp9KhWJl17Nhx6YgYTeEkZWZl1759e3r37l3pMKwK+HBfDVsYanwlM7MKcpKqYQuXVDoCM7MVc5IyM7PCcpIyM7PCcpIyM7PC8tV9Ne6Cu6aUra3tenT2ALNm1iw1naQk7Q9cBrQFromIS+ot3xj4E7BOWuf0iCjfsOEt5NcPPM9lD72wTNnNq2UPL5s2e/nn23RZazW6duqwXPln8co7H7Jwka/UMLPmqdkkJaktcAWwDzALmCBpZERMza12NnBrRPxeUl9gFNCr7MF+Rqfsszmn7LP5soXX/Y7xL7/DTd/ZuSwxXHDXFCcpM2u2Wj4ntRMwPSJeiohPgJuBg+utE8Daaboz8HoZ4zMzq3k125MCegAzc/OzgEH11jkPuF/SCcCawN6lKpJ0HHAc0CJPoiyn9V7+e1naaf9RZ4ha/p/IzFZGLSeppjgSGB4Rv5S0C3CDpG0iYpnjVhFxFXAVZE/mrUCcK23h6l3L0k60BS38uCxtmdmqo5b/tX0N2Cg33zOV5X0buBUgIh4DOgJdyhKdmZnVdJKaAPSR1FvSasARwMh667wK7AUgaSuyJLXKDNvctcPiSodgZrZCNZukImIR8APgPmAa2VV8UyRdIOmgtNr/AN+R9CRwEzAsVqEH4HTrsKjSIZiZrVBNn5NK9zyNqld2Tm56KrBbueMyM7NMzfakzMys+JykzMyssGr6cJ+VVxCMnFz/AsrW0Wn19gzZoltZ2jKz1uMkZWXTvk0bunbqWJa25sxbUJZ2zKx1+XCfmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVlpOUmZkVVk0nKUn7S3pO0nRJpzewzjckTZU0RdKN5Y7RzKyWtat0AJUiqS1wBbAPMAuYIGlkREzNrdMHOAPYLSLeldStMtG2jtU+ebfSIZiZrVAt96R2AqZHxEsR8QlwM3BwvXW+A1wREe8CRMRbZY6xVXVwkjKzgqvlJNUDmJmbn5XK8jYHNpc0TtJ4SfuXLTozM6vdw31N1A7oAwwGegJjJW0bEe/lV5J0HHAcwMYbb1zmEM3MVl213JN6DdgoN98zleXNAkZGxMKIeBl4nixpLSMiroqIgRExsGvXrq0WsJlZranlntQEoI+k3mTJ6QjgqHrr3AkcCVwnqQvZ4b+Xyhlka+s18adlaednH8HYNoOAPcvSnpmtGmo2SUXEIkk/AO4D2gJ/jIgpki4AJkbEyLRsX0lTgcXAaRHxTuWiXkmjL4Z/XFJy0ZrvTluu7JOOXVi4esv2CHsveQViCfNbtFYzW9UpIiodwypl4MCBMXHixEqH0TTndWbKPuW59WvJQz/NktQ37ylLe3PmLeCg/vWvgzGzopI0KSIG1i+v5XNSZmZWcE5SZmZWWDV7TsoqY72X/16ehhZ1ILvmxcyqmZOUlY/atPgFGQ16eyYjJ9e/o8DMWlOn1dszZIuWHT3OSaqGvd59n0qH0GrWW7MjdOpY6TDMasqceQtavE6fk6phr2+4X6VDMDNbIScpK5tPllQ6AjOrNk5SVjYLnaTMrJmcpMzMrLCcpMzMrLCcpMzMrLB8CbqV1emPlqedL3btyM69y9OWmbUeJylrcSOegxufX7bs5tWyn0+XGEO+2+qw/hot1/5L74MWd2DnlqvSzCrEScpa3NAtslder4lZgrr7K63f/umPkj1Yxcyqns9JmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJWdl0W73SEZhZtXGSsrJpycvMzaw2OEmZmVlhOUmZmVlhOUmZmVlhOUmZmVlhOUmZmVlheey+Grb6am1578NPytLWosVLWBJlacrMViFOUjVs6w3XhrW6lKex6e15/6OF5WnLzFYZPtxnZmaF5SRlZmaF5SRlZmaF5SRlZmaFVdNJStL+kp6TNF3S6StY7xBJIWlgOeMzM6t1NZukJLUFrgAOAPoCR0rqW2K9TsBJwL/KG6GZmdVskgJ2AqZHxEsR8QlwM3BwifV+CvwcWFDO4MzMrLaTVA9gZm5+VipbStL2wEYRcfeKKpJ0nKSJkibOmTOn5SM1M6tRtZykVkhSG+BXwP80tm5EXBURAyNiYNeuXVs/ODOzGlHLSeo1YKPcfM9UVqcTsA0wRtIMYGdgpC+eMDMrn1pOUhOAPpJ6S1oNOAIYWbcwIt6PiC4R0SsiegHjgYMiYmJlwjUzqz01m6QiYhHwA+A+YBpwa0RMkXSBpIMqG52ZmUGNDzAbEaOAUfXKzmlg3cHliMnMzD5Vsz0pMzMrPicpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpMzMrrJp+npSV15ofvEqviT9t9XZ+9hGMbTMI2LPV2zKz1uUkZeWx6WA++GhhWX7hei95BWIJ88vQlpm1LicpK4/N9+eZNgNZZ43VWr2pJQ/9lE+WRKu3Y2bLGvX0bA7q36NF6/Q5KVslLVyiSodgVnPunfJmi9fpJGVmZoXlw321rGNnmN/y//k0ZI0PZsMam5atPTOrfk5StazPPmVtbvGMq8vanplVPycpW2VdcNeUSodgZp+Rk5RVtRHPwY3PL1t2c7qAcNrsecut32Wt1ejaqUMZIjNbdc2Z9zFvz/+k5LJep9+9XNlJe/XhlH02X6m2nKSsqg3dInvl9ZoIT78DN31n58oEZVajjrx6PDMu+XKL1umr+8zMrLCcpMzMrLCcpMzMrLCcpMzMrLCcpGyVtH7HxZUOwazm7L/1+i1eZ00nKUn7S3pO0nRJp5dY/kNJUyU9JekhSZtUIk5rvvU7Lql0CGY150vbdm/xOms2SUlqC1wBHAD0BY6U1Lfeak8AAyNiO+AvwP+WN0ozs9pWs0kK2AmYHhEvRcQnwM3AwfkVImJ0RHyYZscDPcsco5lZTavlJNUDmJmbn5XKGvJt4J5WjcjMzJbhESeaQNI3gYHAFxpYfhxwHMDGG29cxsjMzFZttZykXgM2ys33TGXLkLQ3cBbwhYj4uFRFEXEVcBXAwIED/UjYBqy+Wlve+7D0eF8tadHiJbDEF06YrQpqOUlNAPpI6k2WnI4AjsqvIGkA8Adg/4h4q/whrlq23nBtWKtL6zc0vT3zPljY+u2YWaur2XNSEbEI+AFwHzANuDUipki6QNJBabVLgbWA2yRNljSyQuGamdWkWu5JERGjgFH1ys7JTe9d9qDMzGypmu1JmZlZ8TlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYbWrdABWQzp2hvlvtn47ixcSatv67ZhZq3OSsvLps0952pl4HfHJnPK0ZWatyof7zMyssJykzMyssJykzMyssJykzMyssHzhhK2S2rQRc+YtqHQYZjWl0+rtW7xOJylbJXXq0I6D+veodBhm9hn5cJ+ZmRWWk5SZmRWWk5SZmRVWTScpSftLek7SdEmnl1jeQdItafm/JPWqQJhmZjWrZpOUpLbAFcABQF/gSEl96632beDdiNgM+DXw8/JGaWZW22o2SQE7AdMj4qWI+AS4GTi43joHA39K038B9pKkMsZoZlbTajlJ9QBm5uZnpbKS60TEIuB9YL36FUk6TtJESRPnzPHAphW3wbaw7qaVjsLMWoDvk2oBEXEVcBXAwIEDo8Lh2AGXVDoCM2shtdyTeg3YKDffM5WVXEdSO6Az8E5ZojMzs5pOUhOAPpJ6S1oNOAIYWW+dkcAxafpQ4OGIcE/JzKxMavZwX0QskvQD4D6gLfDHiJgi6QJgYkSMBK4FbpA0HfgPWSIzM7MyqdkkBRARo4BR9crOyU0vAA4rd1xmZpap5cN9ZmZWcE5SZmZWWE5SZmZWWE5SZmZWWPIV1S1L0hzglZXcvAvwdguG05qqKVaorngda+uoplihuuJtiVg3iYiu9QudpApE0sSIGFjpOJqimmKF6orXsbaOaooVqive1ozVh/vMzKywnKTMzKywnKSK5apKB9AM1RQrVFe8jrV1VFOsUF3xtlqsPidlZmaF5Z6UmZkVlpOUmZkVlpNUQUjaX9JzkqZLOr3S8eRJ2kjSaElTJU2RdFIqP0/Sa5Imp9eXKh0rgKQZkp5OMU1MZetKekDSC+nn5woQ5xa5fTdZ0lxJJxdpv0r6o6S3JD2TKyu5L5W5PP0OPyVp+wLEeqmkZ1M8d0haJ5X3kvRRbh9fWYBYG/zcJZ2R9utzkvYrQKy35OKcIWlyKm/5/RoRflX4RfaokBeBTYHVgCeBvpWOKxdfd2D7NN0JeB7oC5wHnFrp+ErEOwPoUq/sf4HT0/TpwM8rHWeJ34E3gE2KtF+BPYHtgWca25fAl4B7AAE7A/8qQKz7Au3S9M9zsfbKr1eQ/Vryc09/a08CHYDe6buibSVjrbf8l8A5rbVf3ZMqhp2A6RHxUkR8AtwMHFzhmJaKiNkR8XiangdMA3pUNqpmOxj4U5r+E/DVyoVS0l7AixGxsqOVtIqIGEv2LLW8hvblwcD1kRkPrCOpe1kCpXSsEXF/RCxKs+PJnsBdcQ3s14YcDNwcER9HxMvAdLLvjLJYUaySBHwDuKm12neSKoYewMzc/CwKmgQk9QIGAP9KRT9Ih1L+WIRDaEkA90uaJOm4VLZ+RMxO028A61cmtAYdwbJ/6EXcr3Ua2pdF/z3+FllPr05vSU9I+oekPSoVVD2lPvci79c9gDcj4oVcWYvuVycpazJJawF/BU6OiLnA74HPA/2B2WTd/iLYPSK2Bw4Avi9pz/zCyI5LFObeC0mrAQcBt6Wiou7X5RRtXzZE0lnAImBEKpoNbBwRA4AfAjdKWrtS8SVV87nnHMmy/1y1+H51kiqG14CNcvM9U1lhSGpPlqBGRMTtABHxZkQsjoglwNWU8RDEikTEa+nnW8AdZHG9WXfoKf18q3IRLucA4PGIeBOKu19zGtqXhfw9ljQMOBAYmpIq6dDZO2l6Etl5ns0rFiQr/NyLul/bAV8Hbqkra4396iRVDBOAPpJ6p/+qjwBGVjimpdJx52uBaRHxq1x5/nzD14Bn6m9bbpLWlNSpbprsxPkzZPvzmLTaMcDfKhNhScv8N1rE/VpPQ/tyJHB0uspvZ+D93GHBipC0P/Aj4KCI+DBX3lVS2zS9KdAHeKkyUS6NqaHPfSRwhKQOknqTxfrvcsdXwt7AsxExq66gVfZrua4Q8avRK2i+RHbV3IvAWZWOp15su5Md0nkKmJxeXwJuAJ5O5SOB7gWIdVOyK6GeBKbU7UtgPeAh4AXgQWDdSsea4loTeAfonCsrzH4lS56zgYVk50K+3dC+JLuq74r0O/w0MLAAsU4nO59T93t7ZVr3kPT7MRl4HPhKAWJt8HMHzkr79TnggErHmsqHA8fXW7fF96uHRTIzs8Ly4T4zMyssJykzMyssJykzMyssJykzMyssJykzMyssJymzFiBpHUnfy81vKOkvLVT3eZJOTdMXSNq7hertLumulqirgfrnN2PdBws4/JMVgJOUWctYB1iapCLi9Yg4tKUbiYhzIuLBFqruh2QjGxTBDeT2n1kdJymzlnEJ8Pn0DJ1L03N1noFsWB5Jdyp79tIMST+Q9MM0COd4Seum9T4v6d40MO4jkras34ik4ZIOTdMzJJ0v6XFlz8/aMpWvmQYo/Xdqo6ER9Q8B7k3b3C1puzT9hKRz0vQFkr6Tpk+TNCENgHp+LqZvprYmS/pD3YgDueVdJD0m6cup9zY2rftMbgDSkWQjb5gtw0nKrGWcTvaojf4RcVqJ5duQjXO2I3AR8GFkg3A+Bhyd1rkKOCEidgBOBX7XhHbfjmww3d+nbSAbneDhiNgJGAJcmoaIWioNr/NuRHycih4B9pDUmWwg1t1S+R7AWEn7kg1xsxPZAKg7SNpT0lbA4cBuEdEfWAwMzbWzPnA32fOG7gaOAu5L6/YjG5mAiHgX6CBpvSa8Z6sh7SodgFmNGB3Zs7jmSXof+HsqfxrYLo0wvytwWzZUIpA95K4xt6efk8iSIGTjFR5Udx4L6AhsTPYcsDrdgTm5+UeAE4GXyZLKPpLWAHpHxHOpN7Uv8ERafy2ypLUdsAMwIcW9Op8OONuebPik70fEP1LZBOCPacDiOyNici6Gt4ANyYaJMgOcpMzK5ePc9JLc/BKyv8M2wHuph7Ey9S7m079nAYdExHMr2O4jsuRVZwIwkGww0AeALsB3yJJfXZ0XR8Qf8pVIOgH4U0ScUaKNRWn7/YB/QPYAvfTolC8DwyX9KiKuT+t3THGZLeXDfWYtYx7QaWU3juz5XC9LOgyykecl9VvJ6u4DTkij1yNpQIl1nid71Hdd+5+QDcR6GNkhyEfIDh+OzdX5rdTjQ1IPSd3IekqHpmkkrStpk7pqyR40uKWkH6flm5A9JO9q4Bqyx5LXjbS/ATBjJd+zraKcpMxaQGTP0BmXLga4dCWrGQp8W1LdCO4NXfDQmJ+SHWp7StKUNF8/3g+AFyVtlit+BHgrIj5K0z3TTyLifuBG4DFJTwN/ATpFxFTgbLInIT9F1gvrnmtnMdkFEV9Ml+gPBp6U9ATZuazL0qo7AOPj00e9mwF4FHSzWiXpa8AOEXF2AWK5DBgZEQ9VOhYrFp+TMqtREXFHga6me8YJykpxT8rMzArL56TMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywnKTMzKywKpakJP1R0luSnmnGNoMlhaRjc2X9U9mpjWx7vKSjG1mnv6QvNSGOgZIub2rcjdQ1TNJvW6IuM7NVTSV7UsOB/Vdiu2eAb+TmjwSebGyjiLgyIq5vZLX+QKNJKiImRsSJja1nZmafTcWSVESMBf6zEpu+AnSUtL4kkSW6e+oWSvqOpAmSnpT0V0lrpPLz6npbksZI+rmkf0t6XtIeklYDLgAOlzRZ0uGSdpL0mKQnJD0qaYu0/WBJd+Xq/WOq8yVJJ+Zi+WZqY7KkP0hqm8r/O7X7b2C3ldl/Zma1oFDnpCSdlr7Q67/qH1r7C3AYsCvwOPBxbtntEbFjRPQDpgHfbqC5dhGxE3AycG5EfAKcA9wSEf0j4hbgWWCPiBiQlv2sgbq2BPYDdgLOldRe0lbA4cBuEdEfWAwMldQdOJ8sOe0O9G3i7jEzqzntKh1AXkRcClzahFVvBW4hSw43kSWrOttIuhBYB1gLuK+BOm5PPycBvRpYpzPwJ0l9gADaN7De3RHxMfCxpLeA9YG9gB2ACVmHj9WBt4BBwJiImAMg6RZg8xW8VzOzmlWoJCXpNGBoiUVj8+eAIuINSQuBfYCTWDZJDQe+GhFPShoGDG6gubre12Ia3g8/BUZHxNck9QLGNFJXvj4Bf4qIM/IrSvpqA3WYmVk9hUpSzehJQXb4rVtELE49lTqdgNmS2pMlvNeaEcK8tH2dzrnthzWjHoCHgL9J+nVEvCVp3VT3v4DLJK0HzCU7bNnohR9mZrWokpeg3wQ8BmwhaZakhs4dlRQRj0bEnSUW/YQsEYwjO6fUHKOBvnUXTgD/C1ws6QmamdAjYipwNnC/pKeAB4DuETEbOI/svY8jO29mZmYlKCIqHYOZmVlJhbq6z8zMLM9JyszMCquqk5SkdpLmSLqkXvmZzajjGkkN3quUbtIduJLx7SnpcUmLJB1ab9m9kt6ruym46CSdImmKpGck3SSpY6VjKmVlhtuqpGqK17G2jmqKFcofb1UnKbJL0J8HDtOyl/g1KUlJahsRx6aLHFrDq2RXBd5YYtmlwH+1UrstSlIP4ERgYERsA7QFjqhsVA0azsoNt1Upw6meeIfjWFvDcKonVihzvNWepI4ELiNLBrsApF7V6ukKvRH1N5A0X9IvJT0J7FLXU5LUVtLw1FN4WtIp9bZrk5Zf2NTgImJGRDwFLCmx7CGyS96rRTuy/doOWAN4vcLxlPQZhtuqiGqK17G2jmqKFcofb6Huk2qOdLhpb+D/kY0ucSTwaEScLukHaSiiUtYE/hUR/5PqqSvvD/RIPQUkrZPbph0wAngmIi5q0TdSBSLiNUm/IPtn4CPg/oi4v8JhmVkNqOae1IFko0F8BPwV+GrdAK6NWJzWr+8lYFNJ/ydpf7Ibbev8gRpNUACSPgccDPQGNgTWlPTNykZlZrWgmpPUkcDekmaQjb+3HvDFJmy3ICIW1y+MiHeBfmRDHx0PXJNb/CgwpKgXC5TB3sDLETEnIhaSjXu4ayPbmJl9ZlWZpCStDewBbBwRvSKiF/B9ssQFsDANi9ScOrsAbSLir2QjRWyfW3wtMAq4NZ2TqTWvAjtLWiNdoLIXHinDzMqgKpMU8DXg4TTyeJ2/AV+R1AG4Cniq1IUTK9ADGCNpMvBnYJmBYSPiV8ATwA2SmrTfJO0oaRbZ+Hx/kDQlt+wR4DZgrzQs1H7NiLWsIuJfZI9HeRx4muz35qqKBtWAzzrcVrlVU7yOtXVUU6xQ/ng9LJKZmRVWtfakzMysBjhJmZlZYTlJmZlZYZU9SUnaQNLNkl6UNEnSKEmFfXy6pBnpyr+WrvcMSdMlPVfkiyagusYWk9RR0r8lPZnGGjy/0jE1pJpiheqK17G2jorEGhFle5E9Uv0x4PhcWT9gjzLH0a4Z684AurRw+33JnsbbgewG2ReBtuXcB82Md0+yS/KfqXQsTYhVwFppuj3ZAzB3rnRc1R5rtcXrWFedWMvdkxoCLIyIK+sKIuLJiHgEQNJpkiZIeqouQ0vqJWmapKtT5r5f0upp2YmSpqb1b05l60q6M5WNl7RdKj9P0g2SxpFdRt5V0l9TexMk7ZbWWy+1MUXSNWQfSks7GLg5Ij6OiJeB6cBOrdBOi4gqGlssMvPTbPv0KuQlrNUUK1RXvI61dVQi1nInqW3IRodYjqR9gT5kX9b9gR0k7ZkW9wGuiIitgfeAQ1L56cCAiNiObJQIgPOBJ1LZmcD1uWb6AntHRN3AtL+OiB1TfXUjTJwL/DO1dQew8Wd5ww3oAczMzc9KZdYClA0WPBl4C3ggsvu8CqmaYoXqitexto5yx1qkCyf2Ta8nyG4a3ZIsOUE2JM/kND0J6JWmnwJGKBtHblEq2x24ASAiHgbWUzZCBcDIyMb6g2yon9+mnT0SWFvSWmSHtv6ctr8beLdF36W1uohYHNkAwz2BnSRtU+GQGlRNsUJ1xetYW0e5Yy13kpoC7NDAMgEXR0T/9NosIq5Ny/IjSyzm09HbvwxcQXa+ZIIaH7Log9x0G7JjqXXt9ch1Y1vba8BGufmeqcxaUES8B4ymCp7VU02xQnXF61hbR7liLXeSehjoIOm4ugJJ20naA7gP+FbqzSCph6RuDVWkbGiijSJiNPBjoDOwFvAIMDStMxh4OyLmlqjifuCEXH390+RY4KhUdgDwuZV5o40YCRwhqYOk3mQ9xn+3Qjs1J51rXCdNr072YMxnKxpUA6opVqiueB1r66hErGUdLDUiQtLXgN9I+jGwgOzquZMj4gVJWwGPKXvG03zgm2Q9p1LaAn+W1JmsF3Z5RLwn6Tzgj5KeAj4Ejmlg+xOBK9J67ciS0/Fk57RuUjbO3qNkg6u2qIiYIulWYCrZYcrvR4mR2YtC2Vhdg4EuysYiPDfXyy2a7sCflD22pQ1wa0TcVeGYGlJNsUJ1xetYW0fZY/XYfWZmVlhFunDCzMxsGU5SZmZWWE5SZmZWWFWRpCSNUTbG3eT0+ksL199L0lEtWWcT2qymsfs2kjQ6je4xRdJJlY6pIZK2yP2eTJY0V9LJlY6rlGqKFaorXsfaOioRa1VcOCFpDHBqRExspfoHp/oPbI36S7TXF7iJbHSNDYEHgc2LeoWfpO5A94h4XFInshuqvxoRUysc2gqlK5BeAwZFxCuVjmdFqilWqK54HWvrKFesVdGTKkVSZ0mvpPulkLSmpJmS2kv6vKR7lY2y/oikLdM6wyVdLulRSS9JOjRVdwmwR/rP4BRJWysb6XeysjEA+zQUx0qqtrH7ZkfE42l6HjCN6hjGaS/gxaL/sSfVFCtUV7yOtXWUJdZqSlIjcl3MSyPifWAy8IW0/EDgvohYCFwFnBAROwCnAr/L1dOdbOikA8mSE2RjAD6SRp74Ndn9UpeloT8Gko2t15Kqduw+Sb2AAWSjHxfdEWQ91mpQTbFCdcXrWFtHWWIt6828n9HQEof7bgEOJxua4wjgd8pGrNgVuC3dFAzZIzHq3BkRS4CpktZvoK3HgLMk9QRuj4gXWupNVLO0b/9KdvN1qVE8CkPSasBBwBmVjqUx1RQrVFe8jrV1lDPWaupJlTIS2F/SumRjAj5M9p7ey43J1z8itsptkx8HsORjOCLiRrIP4CNglKQvtnDcVTd2n6T2ZAlqRETcXul4muAA4PGIeLPSgTRBNcUK1RWvY20dZYu1qpNUGhB2AtljN+5Ko/POBV6WdBiAMv0aqWoe0KluRtKmwEsRcTnwN2C7Fg69qsbuU9YlvRaYFhG/qnQ8TXQk1XPYpJpiheqK17G2jrLFWk1X93Un69lANmjs3mnZocBtwOCI+Ecq6w38Pm3TnuwihQskDSdLZn9J682PiLVSL+E+YD1gONnhwf8CFgJvAEdFRIs+9E/SWcC3yMbuOzki7mnJ+luSpN3JBu59GliSis+MiFGVi6phktYkG3Nx03TusrCqKVaorngda+sod6xVkaTMzKw2VfXhPjMzW7U5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWE5SZmZWWGtMElJ2kjSaElTJU2RdFJTK5Y0WFJI+kqu7C5JgxvZ7mRJazS1neaQ9GgT1mlS+5KukdS3heKaIalLS9RlZrYqaawntQj4n4joC+wMfL+ZX8yzgLOaGdPJQKskqYjYtaXaj4hjI2LqZw7KzMwatMIkFRGzI+LxND0PmAb0aEb9TwLvS9qn/gJJe0l6QtLTkv6YHqV+IrAhMFrS6BLbzJB0saTJkiZK2l7SfZJelHR8WmctSQ9JejzVfXBu+/np52BJYyT9RdKzkkakx8wv176k36e2pkg6P1fXGEkD6+qVdJGkJyWNl7R+Ku8q6a+SJqTXbql8PUn3pzqvAdSMfWpmVjsiokkvoBfZI4PXTvOnAZNLvC5PywcDdwF7Av9IZXel8o7ATGDzVH492SPUAWYAXRqIYQbw3TT9a+ApoBPQFXgzlbfLxdgFmM6nTyCen4vtfaAnWaJ+DNi9VPvAuulnW2AMsF2aHwMMTNMBfCVN/y9wdpq+MVfvxsC0NH05cE6a/nLavuR79ssvv/yq5Ve7RnIYkPVOgL+mRDIXICIuBS5tbNuIGCsJSbvnircAXo6I59P8n4DvA79pQjgj08+ngbUi6+HNk/SxpHWAD4CfSdoTWELW81sfeKNePf+OiFnp/U0mS8L/LNHeNyQdR5b8ugN9yZJj3idkCRhgElDXc9wb6Cst7SitnfblnsDXASLibknvNuF9m5nVnEaTlKT2ZAlqRETcnis/DRhaYpOxEXFivbKLgLPJznF9Vh+nn0ty03Xz7VJMXYEdImKhpBlkPbeG6gFYTIl9Iak3cCqwY0S8K2l4A3UtjIgoUVcbYOeIWFCv3gbfnJmZfaqxq/sEXEt2mOpX+WURcWlE9C/xqp+giIj7gc8B26Wi54BekjZL8/8F/CNNzyM7hLeyOgNvpQQ1BNikmdvn21+brGf2fjrPdEAz67ofOKFuRlL/NDkWOCqVHUC2b8zMrJ7Gru7bjSyBfDFdrDBZ0pdWsq2LgI0AUs/iv4HbJD1N1gu6Mq13FXBvqQsnmmgEMDDVezTwbDO3X9p+RDwJPJHquBEY18y6TkyxPCVpKnB8Kj8f2FPSFLLDfq82s14zs5qgT49SmZmZFYtHnDAzs8JykjIzs8Kq6iQlqZ2kOZIuqVd+ZjPqWOHwRvmbdlcivj3TTcWLJB2aK+8v6bF0M+9Tkg5fmfrLSdI6uZufp0napdIxNUTS/pKekzRd0umVjmdFqilWqK54HWvrKHuslb5R67O8yK62Gwe8SDq/lsrnN3H7tk1YZwzppt2ViK8X2RWN1wOH5so3B/qk6Q2B2cA6ld6fjbyXPwHHpunVihov2U3XLwKbpjifBPpWOq5qj7Xa4nWsq06sVd2TAo4ELiO7Om4XgNSrWj1diTii/gZpCKNfSnoS2KWupySpraThkp5JwymdUm+7Nmn5hU0NLiJmRMRTZFcv5sufj4gX0vTrwFtk93YVkqTOZDcgXwsQEZ9ExHsVDaphOwHTI+KliPgEuBk4uJFtKqWaYoXqitexto6yx1q1SUpSR7IRHf4O3ESWsIiI04GPIrtnq9TNxmsC/4qIfhGRH2GiP9AjIraJiG2B63LL2pFd2v5CRJzdwu9jJ7L/SF5syXpbWG9gDnCdsvEWr5G0ZqWDakAPsiG36syieeNNllM1xQrVFa9jbR1lj7VqkxRwIDA6Ij4iGxHjq5LaNmG7xWn9+l4CNpX0f5L2B+bmlv0BeCYiLvqsQedJ6g7cAPx3RCxpbP0KagdsD/w+IgaQ3eBc6OPmZrZqqOYkdSSwdxr2aBKwHvDFJmy3ICIW1y+MiHeBfmTnoI4HrsktfhQYknpvLULS2sDdwFkRMb6l6m0ls4BZEfGvNP8XsqRVRK+RbhpPeqayIqqmWKG64nWsraPssVZlkkpf8HsAG0dEr4joRTZA7ZFplYVpzMHm1NkFaBMRfyUbZzD/JXwtMAq4VVKTBuVtpK3VgDuA6yPiL5+1vtYWEW8AMyVtkYr2Aor6LK0JQB9JvdN+PoJPByUummqKFaorXsfaOsoe62f+wq2QrwEPR0R+kNi/Af8rqQPZ0EZPSXq8gfNSpfQgO+dSl7jPyC+MiF+lCwhukDS0KYfnJO1Ilow+B3xF0vkRsTXwDbILEdaTNCytPiwiJjcx1ko4ARiRfjFfIhvWqnAiYpGkHwD3kV2J9MeImFLhsEqqplihuuJ1rK2jErF6WCQzMyusqjzcZ2ZmtcFJyszMCstJyszMCqvsSUrSBpJulvSipEmSRknavNxxNJWkGenKv5au94w09tVzkvZr6fpbkqQ/SnpL0jOVjqUx1RQrVFe8jrV1VFOsUP54y5qkJInsarcxEfH5iNiB7Cq69cscR0WvalQ2oO0RwNbA/sDvmngjcqUMJ4uzGgynemKF6op3OI61NQynemKFMsdb7p7UEGBhRNQ9hZeIeDIiHgGQdJqkCcpGBj8/lfVSNur21cpGDb9f0upp2YmSpqb1b05l60q6M5WNl7RdKj9P0g2SxpFdRt5V0l9TexMk7ZbWWy+1MUXSNYBaYT8cDNwcER9HxMvAdLIxsQopIsYC/6l0HE1RTbFCdcXrWFtHNcUK5Y+33ElqG7LRIZYjaV+gD9mXdX9gB0l7psV9gCvSPUbvAYek8tOBARGxHcs+mv2JVHYm2QjkdfoCe0dE3cC0v46IHVN9dSNMnAv8M7V1B7DxZ3nDDaimsbrMzCqmSDfz7pteT6T5tciS06vAy7kbXSeRPQID4CmyG0zvBO5MZbuTklhEPJx6RmunZSPTWH+QDU7bNzsCCcDaktYiu8n262n7uyW923Jv0czMmqPcSWoKcGgDywRcHBF/WKZQ6gXkR5ZYDKyepr9MllS+ApwladtG2v8gN90G2DkiFtRrr5EqWkQ1jdVlZlYx5T7c9zDQQdJxdQWStpO0B9kwG99KvRkk9ZDUraGK0vBFG0XEaODHQGey3tcjwNC0zmDg7YiYW6KK+8mG+qmrr3+aHAsclcoOIBvSqKWNBI6Q1EFSb7Ie479boR0zs6pW1iQV2RhMXyMbvfxFSVOAi4E3IuJ+4EbgMUlPk4203WkF1bUF/pzWfQK4PD2I7zyy81lPAZcAxzSw/YnAwHSBxVSWPae1Z4rt62SHG1tUGuvqVrJBWu8Fvl9qZPaikHQT8BiwhaRZkr5d6ZgaUk2xQnXF61hbRzXFCuWP12P3mZlZYXnECTMzKywnKTMzKywnKTMzK6yqSFKSxqQx7ianV4s+zTaNanFUS9bZhDarZuw+WDqG4dNp/08sQDzLjR8m6bA0UsgSSQMrGV9eNcUK1RWvY20dRYq1KpJUMjQi+qdXQ/daraxepMvOy0HVN3ZfnSFp/xfhj2k4y48f9gzZFZljyx7Nig2nemKF6op3OI61NQynILFWU5JahqTOkl5J90shaU1JMyW1l/R5SfcqG2X9EUlbpnWGS7pc0qOSXpJUl+wuAfZIvYRTJG0t6d9p/ilJfVo4/Koau6+ISo0fFhHTIuK5CoXUoGqKFaorXsfaOooUazUlqRG5w32XRsT7wGTgC2n5gcB9EbEQuAo4IY2yfirwu1w93cmGTjqQLDlBNgbgI6mX8Guye6Yui4j+wECysfVaUjWO3RfA/SnxH9fo2mZmLaBIY/c1ZmhE1D8XcgtwODCa7PDZ79KIFbsCt+nTIY465La5MyKWAFMlNfSIkMfIhlnqCdweES+01JuoYrtHxGvKRgF5QNKz6b8tM7NWU009qVJGAvtLWhfYgWzYpTbAe7nzV/0jYqvcNvlxAEsO1BcRNwIHAR8BoyR9sYXjrrqx+yLitfTzLbLR4X140sxaXVUnqYiYD0wge+zGXRGxOI3T97KkwyB70KKkfo1UNY/cEEySNgVeiojLgb8B27Vw6FU1dl8639epbppstPqqeIqomVW5iCj8CxgDPEd2Dmoy8GBu2aFk50u+kCvrTTYm3pNk4+Odk8qHA4fm1puffrYn64U9CZxCdo5qSmrrXmDdVnhPZwEvpvd1QKX3cSOxbpr2zZNpv5xVgJhuAmYDC8nO6X2bbFzIWWS95TfJzlEWYf9VTazVFq9jXfVj9dh9ZmZWWFV9uM/MzFZtTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZYTlJmZlZY/x85fbUZJH544wAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "stride = 16\n",
    "agg_d1, agg_n1 = aggregate(d1, n1, stride)\n",
    "agg_d2, agg_n2 = aggregate(d2, n2, stride)\n",
    "agg_d1, agg_n1, agg_d2, agg_n2 = [list(map(int, await mpc.output(_))) for _ in (agg_d1, agg_n1, agg_d2, agg_n2)]\n",
    "T1_, E1_ = events_from_table(agg_d1, agg_n1)\n",
    "T2_, E2_ = events_from_table(agg_d2, agg_n2)\n",
    "T1_, T2_ = [t * stride for t in T1_], [t * stride for t in T2_]\n",
    "fit_plot(T1_, T2_, E1_, E2_, 'Aggregated Kaplan-Meier survival curves', 'weeks', '1=Maintained', '2=Not maintained')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Picking `stride = 16` achieves a reasonable balance between privacy and utility. To enhance both privacy and utility at the same time, one may look for differentially private randomization techniques, adding a suitable type of noise to the Kaplan-Meier curves. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Secure Logrank Tests\n",
    "\n",
    "The function `logrank_test()` performs a secure logrank test on a secret-shared dataset, similar to function `lifelines.statistics.logrank_test()` used above for a dataset in the clear. The input parameter `secfxp` specifies the secure type to be used for fixed-point arithmetic, and the output is an instance of `lifelines.statistics.StatisticalResult`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.06533921951089436\n"
     ]
    }
   ],
   "source": [
    "print((await logrank_test(secfxp, d1, d2, n1, n2)).p_value)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Relying solely on p-values is in general not a good idea, and this is especially true when handling otherwise (mostly) hidden data. Together with the aggregated curves, however, the p-value may lead to a useful conclusion for a study. \n",
    "\n",
    "The function `logrank_test()` uses one secure fixed-point division per time moment in $1..maxT$. Even though these divisions can all be done in parallel, the total effort is significant when $maxT$ is large. However, \"most of the time\" there is actually no event happening and no divisions need to be performed for these time moments. E.g., in the survival tables for the aml dataset above, there are only 7 time moments with nonzero `d1` entries on the entire timeline $1..161$, and only 9 time moments with nonzero `d2` entries.\n",
    "\n",
    "Therefore, it may be advantageous to first extract the nonzero rows of the survival tables, and then limit the computation of the logrank test to those rows. The extraction needs to be done obliviously, not leaking any information about (the location of) the nonzero entries of the survival tables. To prevent this oblivious extraction step from becoming a bottleneck, however, we will actually exploit the fact that the aggregate curves are revealed anyway. We may simply use `agg_d1` and `agg_d2` to bound the number of events per stride, and extract the nonzero rows obliviously and efficiently for each stride. \n",
    "This is basically what the function `agg_logrank_test()` does:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.06533925483923529\n"
     ]
    }
   ],
   "source": [
    "print((await agg_logrank_test(secfxp, d1, d2, n1, n2, agg_d1, agg_d2, stride)).p_value)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Even for a small dataset like aml, the speedup is already noticeable. For larger datasets, the speedup gets really substantial, as can be noticed for some of the other datasets included with the [kmsurival.py](kmsurvival.py) demo. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Summary\n",
    "We end with two complete runs of the demo on the aml dataset, showing the Chi2 test statistic and p-value for each logrank test. \n",
    "\n",
    "The help message included with the demo shows the command line options:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Showing help message for kmsurvival.py, if available:\n",
      "\n",
      "usage: kmsurvival.py [-h] [-i I] [-s S] [-a A] [--collapse] [--print-tables]\n",
      "                     [--plot-curves]\n",
      "\n",
      "optional arguments:\n",
      "  -h, --help          show this help message and exit\n",
      "  -i I, --dataset I   dataset 0=btrial(default) 1=waltons 2=aml 3=lung 4=dd\n",
      "                      5=stanford_heart_transplants 6=kidney_transplant\n",
      "  -s S, --stride S    interval length for aggregated events\n",
      "  -a A, --accuracy A  number of fractional bits\n",
      "  --collapse          days->weeks->month->years\n",
      "  --print-tables      print survival tables\n",
      "  --plot-curves       plot survival curves\n"
     ]
    }
   ],
   "source": [
    "!python kmsurvival.py -h"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Complete Run:  5 logrank tests + survival curves\n",
    "To show the plots of the survival curves the `main()` function of the demo is called directly from a notebook cell:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using secure fixed-point numbers: SecFxp64:32\n",
      "Dataset: aml, with 1-party split, time 1 to 161 (stride 16) weeks\n",
      "Chi2=3.396389, p=0.065339 for all events in the clear\n",
      "Chi2=3.396389, p=0.065339 for own events in the clear\n",
      "Chi2=2.685357, p=0.101275 for aggregated events in the clear\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZsAAAGcCAYAAAAPnQ1AAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAABFC0lEQVR4nO3debxVZdn/8c+XGQVBwQEBBQ1RVAbFkZznIWmgEi01S5+sHBos03459PRkj1aPlZUzzmY4REppgyQpKg6IoDkhKmaIBAgKCnj9/lj3gc1hH87eh73O3hu+79frvM7ea7jXtdc5e1973Wut61ZEYGZmlqc21Q7AzMzWfU42ZmaWOycbMzPLnZONmZnlzsnGzMxy52RjZma5c7Kx9YKkkyT9o9pxAEiaLmn/CrQzU9LBLVz3eEn3t3DdfpJCUrv0/I+STixhvX0kPb+G+WMk/XdLYlpTfFYbnGzWQelDaLGkRZJmpzdxlxa2NUHSl9Yiln6SHpD0nqR/tvTDMW+SvpjiW5j22XhJXfPYVkTsGBET8mi7jBhujohDK9TWERFxfQnLTYyIgZXYZr2ppS871eJks+76WER0AXYBhgPfK2dlZSrx/3Er8BTQAzgPGCtp0wq0WzGS9gP+BxgdEV2BHYDftrAtf5s2K8LJZh0XEW8AfwR2krSxpHskzZE0Lz3u07BsOor5oaSHgPeAG4F9gF+mo6RfSrpc0k8KtyFpnKSvN962pO3Ikt35EbE4Iu4AngE+VUrsko6S9JSkdyS9LumCgnkNXSVfSPPmSfqypN0kTZU0X9IvS9xNuwGTIuKptM/+ExHXR8TCgv2y4uiu8bfUFMdXJb0IvCjp15IubfRafi/pG+nxTEkHS9oyHYFuUrDcMElvS2ovaVtJf5M0N027WVL3El/TGjXxGr4s6cW07y6XpDSvraRLUwwzgKMatTVB0pckdUzr7lQwb9P0GjeTtL+kWY1e65PpaPK3QKem4iuI8SPpcZP/GyW89i0l3ZHeB69IOqNgepN/j/T8ZEnPpf+3+yRt3dw+lLQD8Btgr/Q+mp+WP1LSs+n1vyHpW6W+hnrkZLOOk9QXOJLs6KINcB2wNbAVsBho/IH8eeBUoCtwEjAR+FpEdImIrwHXA6OVjnok9QQOBm4psvkdgRkNH9rJ02k6kj7a8MZrwrvACUB3sg+40yR9vNEyewADgM8C/0d29HRw2sZnlB21NOdR4DBJF0oaIaljCes09vEUyyCyo7nPFnxYbwwcCtxWuEJE/AuYxKrJ9zhgbEQsBQT8CNiS7GirL3BBC2Ir1dFkiXcw8BngsDT9lDRvGNlR8qhiK0fE+8CdwOiCyZ8B/h4RbxUuK6kDcDfZF5pNgN9R4peQpJT/jdWk/9s/kP0f9gYOAs6SdFhzfw9JI4FzgU8Cm5K9N25ttInV9mFEPAd8mewLTZeI6J6WvQb4r3Q0vRPwtzJef91xsll33Z0+yP8B/B34n4iYGxF3RMR7KQH8EGj8YTwmIqZHxLL0gbeKiHgMWED2JgU4FpgQEbOLxNAlLVtoAVkiIyL+UfDGW01ETIiIZyLiw4iYSvbGbhzvDyJiSUTcT/YBdGtEvJWO6CaSfUCuUURMJPsA2QW4F5gr6aeS2ja3boEfpSOixWm7QXZUCNmH86T0YdbYLaQP55Scjk3TiIiXIuLPEfF+RMwBflrk9VfSxRExPyJeAx4AhqbpnwH+LyJej4j/kCXAptxC9hoaHEfxLyJ7Au1Tu0sjYiwwudRAS/zfKGY3YNOIuCgiPoiIGcBVBTE3+fcgSxg/iojnImIZWdfr0MKjG5reh8UsBQZJ2igi5kXEk6W89nrlZLPu+nhEdI+IrSPiKxGxWNIGkq6Q9Kqkd4AHge6NPlRfL6Ht64HPpcefI/t2WswiYKNG0zYCFhZZdjWS9lB2ccEcSQvI3uw9Gy1WmOQWF3le0oUREfHHiPgY2bfskWRHdeVcGLFiv0VW3fY2Vn7DPw64uYn17iDrXukF7At8SJaskLS5pNtSF8s7wE2s/vpXI2mr1F2zSNKiMl7Dvwsev8fKfbclq/5fvLqGNh4ANkh/u35kH7Z3FVluS+CNWLUS8JraXUWJ/xvFbA1smbq55qcvZOcCm6f5Tf490rqXFaz3H7Kjz94F7Te1D4v5FFmvw6uS/i5prxLir1tONuuXbwIDgT0iYiOyNxNkb5gGjcuAFysLfhMwUtIQsu6du5vY3nRgG616VdeQNL0UtwDjgL4R0Y2s31trXmXtpG/KfyXr0mg49/AusEHBYlsUW7XR81uBUelb7x5kH2LFtjcPuJ+sG/A44LaCD+D/Se3unP5en6OE1x8Rr6Xumi7pIpG19SZZF16Drdaw7eXA7WSJdjRwT6Nu1MI2ezd0NRZpd5V9LqnxPm/p/8brwCvpi1jDT9eIODLFv6a/x+tk3V6F63aOiIdL2O5q76OImBwRI4HNyN5Dt5fQTt1yslm/dCX7tj8/nQQ9v4R1ZgPbFE6IiFlkXR43AnekrqPVRMQLwBTgfEmdJH2CrC+76AdvE/H+JyKWSNqd7M1fcZJGSjpW2QUUStvaD3gkLTIF+GQ6MvwI8MXm2kwXG7wNXA3cFxHz17D4LWTnH0axapdTV7KjwwWSegNnl/fKKuZ24AxJfdL5p3OaWf4Wsg/r4ynehQbZuZFlqd32kj4J7F4w/2lgR0lDJXVi9XNVLf3feAxYKOk7kjoru/hhJ0m7NYq/2N/jN8B3JTWcc+wm6dMlbnc20Cedq0JSB2X3OnVL3dXvkB1FrbOcbNYv/wd0JvsQfAT4UwnrXEb2DX2epJ8XTL8e2Jmmu9AaHEt2UnkecDEwKp1/aLjJb03dPF8BLpK0EPg++X3zm0d2EvxFsjf9TcAlEdHQ9fUz4AOyD4zrabpLrLFbaPriiULjyC5y+HdEPF0w/UKy80gLyM4l3VnidivtKuA+sgTwZHNxRMSjZEcmW5JdCVlsmQ/IzpOdRNYd9dnCdtMXlYuAv5D9XRrfo9Ki/4105HU0WffeK6z8QtCtYLGif4+IuAv4MXBb6tacBhxRynbJjpSnA/+W9Haa9nlgZmrry2TJeZ2l8OBp1gKS9iX7UN46/E9kZs3wkY2VTdk9B2cCVzvRmFkpnGysLOkGtflAL7JuOTOzZrkbzczMcucjGzMzy52TjZmZ5c4Vaovo2bNn9OvXr9phmJnVlSeeeOLtiCha1d3Jpoh+/frx+OOPVzsMM7O6IqnJkkPuRjMzs9w52ZiZWe6cbMzMLHc+Z2NmrW7p0qXMmjWLJUuWVDsUa4FOnTrRp08f2rdvX/I6dZ1sJF1LVlTvrYjYqch8kRWSPJJsbImT1vUBiszqwaxZs+jatSv9+vVj1VEGrNZFBHPnzmXWrFn079+/5PXqvRttDHD4GuYfQVa9dQDZUMe/boWYzKwZS5YsoUePHk40dUgSPXr0KPuotK6TTUQ8SFaevCkjgRsi8wjZqJS9Wic6M1sTJ5r61ZK/XV13o5WgN6sOZzsrTXszj4098qtT6Dr/uZKWfajzAfx1gyOLzhs5tDfH7dHkYIhmVgGSOP7447npppsAWLZsGb169WKPPfbgnnvuaXK9xx9/nBtuuIGf//znTS4zf/58brnlFr7yla80G8fee+/Nww+XMtjnms2cOZOjjz6aadOmrXVbeVjXk03JJJ1K1tXGVlu1/IP+w2XLml2m/4ev8uGyZdzV5pDV5r069z3eWbLUycYsZxtuuCHTpk1j8eLFdO7cmT//+c/07t272fWGDx/O8OHD17jM/Pnz+dWvflVSsqlEoqkH63qyeYNVx07vk6atJiKuBK4EGD58eItKYe/5latKW/C6o9h5yXz+eNq+q8367BWTeGfJ0pZs3szKdOSRR3LvvfcyatQobr31VkaPHs3EiRMBeOyxxzjzzDNZsmQJnTt35rrrrmPgwIFMmDCBSy+9lHvuuYcLLriA1157jRkzZvDaa69x1llnccYZZ3DOOefw8ssvM3ToUA455BDOP/98Ro4cybx581i6dCn//d//zciRIwHo0qULixYtYsKECVxwwQX07NmTadOmseuuu3LTTTchiSeeeIJvfOMbLFq0iJ49ezJmzBh69erFE088wcknnwzAoYceWrX9WIp1PdmMA74m6TZgD2BBROTShVZJyz8Mxk0pmhMB6Nq5PQcM3KwVIzLLz4V/mM6z/3qnom0O2nIjzv/Yjs0ud+yxx3LRRRdx9NFHM3XqVE4++eQVyWb77bdn4sSJtGvXjr/85S+ce+653HHHHau18c9//pMHHniAhQsXMnDgQE477TQuvvhipk2bxpQpU4Csi+6uu+5io4024u2332bPPffkmGOOWe3cx1NPPcX06dPZcsstGTFiBA899BB77LEHp59+Or///e/ZdNNN+e1vf8t5553Htddeyxe+8AV++ctfsu+++3L22Wev/Y7LUV0nG0m3AvsDPSXNAs4H2gNExG+A8WSXPb9EdunzF6oTaXk+/DDYtGunJufPWeh7E8wqYfDgwcycOZNbb72VI49c9RzqggULOPHEE3nxxReRxNKlxXscjjrqKDp27EjHjh3ZbLPNmD179mrLRATnnnsuDz74IG3atOGNN95g9uzZbLHFFqsst/vuu9OnTx8Ahg4dysyZM+nevTvTpk3jkEOybvfly5fTq1cv5s+fz/z589l336yH5POf/zx//OMf13qf5KWuk01EjG5mfgBfbaVwzKwFSjkCydMxxxzDt771LSZMmMDcuXNXTP9//+//ccABB3DXXXcxc+ZM9t9//6Lrd+zYccXjtm3bsqzIedubb76ZOXPm8MQTT9C+fXv69etX9NLhYm1FBDvuuCOTJk1aZdn58+eX+Uqrq64vfTYzW1snn3wy559/PjvvvPMq0xcsWLDigoExY8aU1WbXrl1ZuHDhKm1tttlmtG/fngceeIBXX22yOPJqBg4cyJw5c1Ykm6VLlzJ9+nS6d+9O9+7d+cc//gFkCa2WOdmY2XqtT58+nHHGGatN//a3v813v/tdhg0bVvRoZU169OjBiBEj2GmnnTj77LM5/vjjefzxx9l555254YYb2H777Utuq0OHDowdO5bvfOc7DBkyhKFDh664gu26667jq1/9KkOHDiXryKldqvUAq2H48OGR63g21x0FS+bDaQ+tNuuzV0xi7qL3+cHHd159vWTOwiUcM7T5SzTNatVzzz3HDjvsUO0wbC0U+xtKeiIiil4XXtfnbOpaLIdnxq4+ffEGwAatHo6ZWZ6cbKrlww+hy+ZFps9rdtUly5av8dJo8OXRZlZbnGzqUN+NN2x2GV8ebWa1xBcImJlZ7pxszMwsd042ZmaWOycbM1svSeKb3/zmiueXXnopF1xwwRrXufvuu3n22WdziedLX/pSs22Xuv3f/OY33HDDDRWJ66STTmLs2CJXzpbJycbM1ksdO3bkzjvv5O233y55nTyTzdVXX82gQYMqsv0vf/nLnHDCCZUKrSKcbMxsvdSuXTtOPfVUfvazn602b+bMmRx44IEMHjyYgw46iNdee42HH36YcePGcfbZZzN06FBefvnlVdY56aSTOO2009hzzz3ZZpttmDBhAieffDI77LADJ5100orlTjvtNIYPH86OO+7I+eefv2L6/vvvT8PN5F26dOG8885jyJAh7LnnnsyePbvo9q+66ip22203hgwZwqc+9Snee+89AC644AIuvfTSFe1+5zvfYffdd2e77bZbUdV6+fLlnH322ey2224MHjyYK664AsiKhn7ta19j4MCBHHzwwbz11luV2d8VacUq6o35i7nonulNzh+xbU8O2qHIPTpm9eiP58C/n6lsm1vsDEdc3OxiX/3qVxk8eDDf/va3V5l++umnc+KJJ3LiiSdy7bXXcsYZZ3D33XdzzDHHcPTRRzNq1Kii7c2bN49JkyYxbtw4jjnmGB566CGuvvpqdtttN6ZMmcLQoUP54Q9/yCabbMLy5cs56KCDmDp1KoMHD16lnXfffZc999yTH/7wh3z729/mqquu4nvf+95q2+/evTunnHIKAN/73ve45pprOP3001eLa9myZTz22GOMHz+eCy+8kL/85S9cc801dOvWjcmTJ/P+++8zYsQIDj30UJ566imef/55nn32WWbPns2gQYNWjJmzNnxkU2NG9ltO7+6dm5z/6tz3eOjl0g/7zaxpG220ESeccMJqQzxPmjSJ4447DshK9zcUu2zOxz72MSSx8847s/nmm7PzzjvTpk0bdtxxR2bOnAnA7bffzi677MKwYcOYPn160W6xDh06cPTRRwOw6667rli3sWnTprHPPvuw8847c/PNNzN9evEvqZ/85CdXa+v+++/nhhtuYOjQoeyxxx7MnTuXF198kQcffJDRo0fTtm1bttxySw488MCSXntzfGRTY44bsJwu2wxocjybNR3xmNWlEo5A8nTWWWexyy678IUvrP1wVw1DBLRp02aV4QLatGnDsmXLeOWVV7j00kuZPHkyG2+8MSeddFLRoQbat2+/YmC1poYtgKzr7u6772bIkCGMGTOGCRMmrDGuwrYigl/84hccdthhqyw7fvz48l50iXxks45qKGnT3M8Dz1emP9asXm2yySZ85jOf4Zprrlkxbe+99+a2224DstL9++yzD7D60AHleuedd9hwww3p1q0bs2fPLnuws8bbX7hwIb169WLp0qVlDzFw2GGH8etf/3rFoHAvvPAC7777Lvvuuy+//e1vWb58OW+++SYPPPBAWe02xclmHdV34w3ZtGunZn8WLi4++qDZ+uSb3/zmKlel/eIXv+C6665j8ODB3HjjjVx22WVANoz0JZdcwrBhw1a7QKAUQ4YMYdiwYWy//fYcd9xxjBgxoqz1G2//Bz/4AXvssQcjRowoa9gCyC61HjRoELvssgs77bQT//Vf/8WyZcv4xCc+wYABAxg0aBAnnHACe+21V1ntNsVDDBTRKkMMvDsHjvrJ6vMWzWbc8r2a7Ub7/tGVGd3QwxVYNXiIgfpX7hADPrIxM7PcOdmYmVnunGzMzCx3TjZmVhU+X1y/WvK38302tWbpErZ8YzwbdW6/2qxlHTYCerZ+TGYV1qlTJ+bOnUuPHj1W3E9i9SEimDt3Lp06Fb+IqSlONrVm4615f+7bLO3cYbVZ7RfPwcnG1gV9+vRh1qxZzJkzp9qhWAt06tSJPn36lLWOk42Ztbr27dvTv3//aodhrcjnbOrQnIXvVzsEM7OyONnUobcXfVDtEMzMyuJutGpZ8Br86ZzVp2+zP7QpegNuLhpqqLVE187tOWDgZhWOyMzWRU421bDzqKxcTWP/eSX7/ZHWSzZ9N96wxevOWbh6tVozs2KcbKph+BegY1fo0mgAtGJHOmZm6wAnmzrV3Lg2Hs3TzGpJ3ScbSYcDlwFtgasj4uJG87cCrge6p2XOiYh8RgeqsJufh1teKJyyKZCNZfHcm6uPqdGzSwc27dqRV+e+B7ztZGNmNaOuk42ktsDlwCHALGCypHERUTjO6veA2yPi15IGAeOBfq0ebAscPzD7adB+8Rzm9v8Yo696hFtP2bPJ9Tyap5nVmnq/9Hl34KWImBERHwC3ASMbLRPARulxN+BfrRifmZlR50c2QG/g9YLns4A9Gi1zAXC/pNOBDYGDizUk6VTgVICtttqq4oFWQpvlS+jxyh+ATdPv4tov7gawxmUaLOuwEQt671epEM3Miqr3ZFOK0cCYiPiJpL2AGyXtFBEfFi4UEVcCV0I2UmcV4mzW+136rni8tPOmTS4XbZtfpkFWb83MLF/13o32BtC34HmfNK3QF4HbASJiEtAJV7M0M2tV9Z5sJgMDJPWX1AE4FhjXaJnXgIMAJO1Almzq+uv8cdtVOwIzs/LUdbKJiGXA14D7gOfIrjqbLukiScekxb4JnCLpaeBW4KSo81GbCq9QMzOrB3V/zibdMzO+0bTvFzx+FhjR2nGZmdlKdX1kY2Zm9cHJxszMcudkY2Zmuav7czZWPU2NheNxbsysMScba7GmxsLxODdm1pi70czMLHdONmZmljsnGzMzy53P2aznVlaSXp0rQptZpTjZrOcKK0k35orQZlYpTjbrqBkL4JyHm19uv95wxNb5x2Nm6zcnm3XQfr1LW27Gguy3k42Z5c3JZh10xNalJZBSjnzMzCrBV6OZmVnunGzMzCx3TjZmZpY7JxszM8udk42ZmeXOycbMzHLnZGNmZrnzfTa15j+vsNMzP6Jd29W/ByzYYm/m9TmoCkGZma0dJ5tass3+2e/FS1eb1WnhqwBONmZWl5xsasl2h8N2hzPtpbfpvkGHVWb1e/wHVQqqfE0NF11pHn7arH442VjFNTVcdKV5+Gmz+uFkY01qaqwbj3NjZuVysrEmNYx1c/PzcPzAldM9zo2ZlcuXPluzbnmh2hGYWb1zsjEzs9w52ZiZWe6cbMzMLHd1f4GApMOBy4C2wNURcXGRZT4DXAAE8HREHNeqQdawGQtKG7GzcBkt78bS6dMZsW1PDtph8/yCM7N1Rl0nG0ltgcuBQ4BZwGRJ4yLi2YJlBgDfBUZExDxJ6+xdgJu+PJY5244qefn9eq/6fPZ78Nbi4ss+M7fwWQeYv5A5C993sjGzktR1sgF2B16KiBkAkm4DRgLPFixzCnB5RMwDiIi3Wj3KVrLZjDvLSjZHbJ39NOeoP8C9H1v5vP3iOZw5fZsWRGhm66t6P2fTG3i94PmsNK3QdsB2kh6S9EjqdjMzs1ZU70c2pWgHDAD2B/oAD0raOSLmFy4k6VTgVICtttqqlUM0M1u31fuRzRtA34LnfdK0QrOAcRGxNCJeAV4gSz6riIgrI2J4RAzfdNNNcwvYzGx9VO9HNpOBAZL6kyWZY4HGV5rdDYwGrpPUk6xbbUZrBlkpnRa+2mz153KrQ3uMHDNrDXWdbCJimaSvAfeRXfp8bURMl3QR8HhEjEvzDpX0LLAcODsi5jbdam1asMXeKx63XzyHDkveLrrchvOeW23aB516srTz6kdrpY6Rc9x25URqZra6uk42ABExHhjfaNr3Cx4H8I30U7fm9Tmo2aSw45+PY/oht5TcZqlHQYVFOM3MWqLez9mYmVkdcLIxM7Pc1X03mrW+NsuX0H5xdtqrcHA1D6pmZk1xsrGyvd+lL9E2e1x44YEHVTOzprgbbR3y1jafrHYIZmZFOdmsQ8qpi2Zm1pqcbMzMLHdONmZmljsnGzMzy52TjZmZ5c6XPlvdWrJsOeOmNC7ybWZro2vn9hwwsPIDGjvZWN3qu/GG1Q7BbJ0zZ+GSXNp1N5qZmeXORzY1aIOO7Zj/3gdNzv9g+Yds1rVTK0ZkZrZ2nGxq0LC+3dc4/6GXio9lY2ZWq5xsrMVmLIBzHl75XMu7sXT6dEZs25ODdti8eoGZWc1xslnPrWmo6TUNGb1f7+LtvTr3PeBtJxszW4WTzXqscKjpxpobMvqIrbOfQu0XL+DM6T0qFp+ZrTucbKqlUzdYNHv16UuXwMZbrz49B2saarrUIaPNzErhZFMtAw4pPv2Zsa0bh5lZK/B9NmZmljsnGzMzy5270axi2ixfQvvFcwHo8cof1qqtZR02YkHv/SoRlpnVACcbq5j3u/Ql2maPl3bedK3aar94TgUiMrNa4W40MzPLnY9s6lBztdPA9dPMrLY42dSh5mqngeunmVltcTeamZnlzsnGzMxy52RjZma5c7IxM7Pc1X2ykXS4pOclvSTpnDUs9ylJIWl4a8ZnZmZ1nmwktQUuB44ABgGjJQ0qslxX4Ezg0daN0MzMoM6TDbA78FJEzIiID4DbgJFFlvsB8GNgSWsGZ2ZmmXpPNr2B1wuez0rTVpC0C9A3Iu5dU0OSTpX0uKTH58xxqRQzs0pap2/qlNQG+ClwUnPLRsSVwJUAw4cPj3wjqw9NDRm9puGi8zL2idcZtWvfVt2mmVVOvR/ZvAEUfgL1SdMadAV2AiZImgnsCYzzRQLNW7DF3izpuvqIoZ0Wvkq3fz/c6vHc8eQbzS9kZjWr3o9sJgMDJPUnSzLHAsc1zIyIBUDPhueSJgDfiojHWznOutPUkNEeLtrMWqKuj2wiYhnwNeA+4Dng9oiYLukiScdUNzozM2tQ70c2RMR4YHyjad9vYtn9WyMmMzNbVd0nG6s9MxbAOWt5WkfLu7F0+vRVpl10z6rPR2zbk4N22HztNmRmrcLJxipqv97NL7Mms9+DtxYDdID5C1eZ99ybqz+f994HvkrNrA442VhFHbF19rO22i+ew9z+H1vxfPRVj3DrKXuueN5wlONEY1Yf6voCATMzqw8+sllHlTJ0NNTu8NFtli+hxyt/KJiyadViMbO152Szjipl6Gio3eGj3+/i7jGzdYm70czMLHc+sqk1nbrBotnF5y1dAhtX4Ox7Hfpc/3erHYKZrQUnm1oz4JCm5z0ztvXiqDEnbPsec6sdhJm1mLvRzMwsd042ZmaWOycbMzPLnZONmZnlzsnGzMxy52RjZma586XPVrZOC18tOmLngi32Ljq6p5mZk42VZcEWexed3mnhqwBONmZWlJONlWVen4OKJpRiRzpmZg18zsbMzHLnZGNmZrlzN1o9WVORzhba4N03YYNtKtqmmVljTjb1ZE1FOluow7+uLWmQtWJqdeA1M6s9TjbruR233Ai69GzRurU68JqZ1R6fszEzs9w52ZiZWe6cbMzMLHdONmZmljsnGzMzy52TjZmZ5c7JxszMclf3yUbS4ZKel/SSpHOKzP+GpGclTZX0V0lbVyNOM7P1WV0nG0ltgcuBI4BBwGhJgxot9hQwPCIGA2OB/23dKM3M6sf4Z97Mpd26TjbA7sBLETEjIj4AbgNGFi4QEQ9ExHvp6SNAn1aO0cysbvxpemXrLzao92TTG3i94PmsNK0pXwT+mGtEZma2mvWmNpqkzwHDgf2amH8qcCrAVltt1YqRrTuaGi66JTzEtNm6pd6TzRtA34LnfdK0VUg6GDgP2C8i3i/WUERcCVwJMHz48Kh8qDVqTcMWLF0CG5d2PUVTw0W3KCQPMW22zqn3ZDMZGCCpP1mSORY4rnABScOAK4DDI+Kt1g+xxq1p2IJnxpbcTFPDRbdEqUdHr859j4vumV6RbZpZvuo62UTEMklfA+4D2gLXRsR0SRcBj0fEOOASoAvwO0kAr0XEMVUL2ipixLY9AQ9xYNZScxa+z9uLio9l1e+ce1ebduZBA/j6Idu1eHt1nWwAImI8ML7RtO8XPD641YOy3B20w+YctMPm1Q7DbJ0z+qpHmHnxURVvt96vRjMzszpQ90c2Vj0bdGxXdEhpDxdtZo052ViLDevbveh0DxdtZo25G83MzHLnZGNmZiscvmM+F9442ZiZ2QpH7twrl3adbMzMLHdONmZmljsnGzMzy50vfbamNVWks4wCnWZm4GRja9JUkc4yCnSamYG70czMrBU42ZiZWe6cbMzMLHc+Z2M1qfEQ04t6DGZu/49VMSIzWxtONlZzGg8x3Wnhq2j50ipFY2aV4GRjFdfU0AOlmr/JPrDJPiue7/TMj2jz4VLmLFyyynJLli2n78Ybtng7ZtZ6nGys4poaeqDFXmoPy+GYob1XmTxuyhuV3Y6Z5cYXCJiZWe6cbMzMLHfuRrPyuYyNmZXJycbK5zI2ZlYmd6OZmVnunGzMzCx3TjZmZpY7JxszM8udk42ZmeXOycbMzHLnS5+tcpq6/6YlfM+O2TrFycYqp6n7b1rC9+yYrVPcjWZmZrmr+2Qj6XBJz0t6SdI5ReZ3lPTbNP9RSf2qEKaZ2XqtrrvRJLUFLgcOAWYBkyWNi4hnCxb7IjAvIj4i6Vjgx8BnWz9aq7SunduvNsaNma2drp3b59JuXScbYHfgpYiYASDpNmAkUJhsRgIXpMdjgV9KUkREawZqZSq82GD5Umiz+kH4AQM3a+WgzKyl6j3Z9AZeL3g+C9ijqWUiYpmkBUAP4O3ChSSdCpwKsNVWW+UVr5Wq8GKDWY/DOx4ozaye1XuyqZiIuBK4EmD48OE+6qklR1xc7QjMbC3V+wUCbwB9C573SdOKLiOpHdANmNsq0ZmZGVD/yWYyMEBSf0kdgGOBcY2WGQecmB6PAv7m8zVmZq2rrrvR0jmYrwH3AW2BayNiuqSLgMcjYhxwDXCjpJeA/5AlJDMza0V1nWwAImI8ML7RtO8XPF4CfLq14zIzs5XqvRvNzMzqgJONmZnlzsnGzMxy52RjZma5k68CXp2kOcCrLVy9J42qE9QIx1Uex1W6WowJHFe5KhHX1hGxabEZTjYVJunxiBhe7Tgac1zlcVylq8WYwHGVK++43I1mZma5c7IxM7PcOdlU3pXVDqAJjqs8jqt0tRgTOK5y5RqXz9mYmVnufGRjZma5c7IxM7PcOdlUkKTDJT0v6SVJ51Qxjr6SHpD0rKTpks5M0zeR9GdJL6bfG1chtraSnpJ0T3reX9KjaZ/9Ng0V0doxdZc0VtI/JT0naa8a2VdfT3+/aZJuldSpGvtL0rWS3pI0rWBa0f2jzM9TfFMl7dLKcV2S/o5TJd0lqXvBvO+muJ6XdFhrxlUw75uSQlLP9Lyq+ytNPz3ts+mS/rdgemX3V0T4pwI/ZEMcvAxsA3QAngYGVSmWXsAu6XFX4AVgEPC/wDlp+jnAj6sQ2zeAW4B70vPbgWPT498Ap1UhpuuBL6XHHYDu1d5XZMOZvwJ0LthPJ1VjfwH7ArsA0wqmFd0/wJHAHwEBewKPtnJchwLt0uMfF8Q1KL0nOwL903u1bWvFlab3JRsO5VWgZ43srwOAvwAd0/PN8tpfPrKpnN2BlyJiRkR8ANwGjKxGIBHxZkQ8mR4vBJ4j+/AaSfbBSvr98daMS1If4Cjg6vRcwIHA2CrG1I3sTXgNQER8EBHzqfK+StoBndMIsxsAb1KF/RURD5KNBVWoqf0zErghMo8A3SX1aq24IuL+iFiWnj5CNnpvQ1y3RcT7EfEK8BLZe7ZV4kp+BnwbKLwqq6r7CzgNuDgi3k/LvFUQV0X3l5NN5fQGXi94PitNqypJ/YBhwKPA5hHxZpr1b2DzVg7n/8jebB+m5z2A+QUfDtXYZ/2BOcB1qXvvakkbUuV9FRFvAJcCr5ElmQXAE1R/fzVoav/U0vvgZLKjBqhyXJJGAm9ExNONZlV7f20H7JO6Zv8uabe84nKyWYdJ6gLcAZwVEe8UzovsWLnVrnuXdDTwVkQ80VrbLFE7sq6FX0fEMOBdsm6hFVp7XwGkcyAjyZLhlsCGwOGtGUOpqrF/miPpPGAZcHMNxLIBcC7w/eaWrYJ2wCZkXXhnA7enHoeKc7KpnDfI+mQb9EnTqkJSe7JEc3NE3Jkmz244RE+/32pq/RyMAI6RNJOsi/FA4DKyboOGEWOrsc9mAbMi4tH0fCxZ8qnmvgI4GHglIuZExFLgTrJ9WO391aCp/VP194Gkk4CjgeNTIqx2XNuSfWl4Ov3/9wGelLRFleOC7P//ztSN9xhZr0PPPOJysqmcycCAdLVQB+BYYFw1AknfTK4BnouInxbMGgecmB6fCPy+tWKKiO9GRJ+I6Ee2b/4WEccDDwCjqhFTiuvfwOuSBqZJBwHPUsV9lbwG7Clpg/T3bIirqvurQFP7ZxxwQrrKak9gQUF3W+4kHU7WVXtMRLzXKN5jJXWU1B8YADzWGjFFxDMRsVlE9Ev//7PILuD5N1XeX8DdZBcJIGk7sgtk3iaP/ZXXlQ/r4w/ZlSUvkF25cV4V4/goWbfGVGBK+jmS7BzJX4EXya5A2aRK8e3PyqvRtkn/xC8BvyNdFdPK8QwFHk/7625g41rYV8CFwD+BacCNZFcGtfr+Am4lO2+0lOyD8otN7R+yq6ouT++BZ4DhrRzXS2TnGhr+739TsPx5Ka7ngSNaM65G82ey8mq0au+vDsBN6X/sSeDAvPaXy9WYmVnu3I1mZma5c7IxM7PcOdmYmVnunGzMzCx3TjZmZpY7JxuzRFn1568UPN9S0tg1rVNG2xdI+lZ6fJGkgyvUbi+lCtp5kLSojGX/oipUx7b64GRjtlJ3YEWyiYh/RcSophdvmYj4fkT8pULNfQO4qkJtra0bKdh/ZoWcbMxWuhjYVtKUNC5Kv4axPySdJOnuNHbLTElfk/SNVLzzEUmbpOW2lfQnSU9Imihp+8YbkTRG0qj0eKakCyU9KemZhuUlbZjGH3ksbaOpCuKfAv6U1rlX0uD0+ClJ30+PL5J0Snp8tqTJaeyUCwti+lza1hRJV0hq2yjmnpImSToqHU09mJadJmmftNg4YHQL972t45xszFY6B3g5IoZGxNlF5u8EfBLYDfgh8F5kxTsnASekZa4ETo+IXYFvAb8qYbtvR8QuwK/TOpDdvf23iNidrJzIJaka9QqpjMi8SOXhgYlkFXy7kRWhHJGm7wM8KOlQsrIju5NVTdhV0r6SdgA+C4yIiKHAcuD4gu1sDtwLfD8i7gWOA+5Lyw4hu1OfiJgHdJTUo4TXbOuZds0vYmbJA5GND7RQ0gLgD2n6M8DgVGV7b+B3BYVzO5bQbkOh1CfIkhlkg4Ad03CeB+gEbEU2NlGDXmTDIzSYCJxBNujavcAhqeJw/4h4Ph3dHAo8lZbvQpZ8BgO7ApNT3J1ZWVizPVlZmq9GxN/TtMnAtanY690RMaUghrfIqlTPLeF123rEycasdO8XPP6w4PmHZO+lNmTjzQxtYbvLWfmeFPCpiHh+DestJktCDSYDw4EZwJ/JqveeQpbEGtr8UURcUdiIpNOB6yPiu0W2sSytfxjwd8gG4ZK0L9lAeGMk/TQibkjLd0pxma3C3WhmKy0kG0a7RSIbM+gVSZ+GFePLD2lhc/cBp6eKz0gaVmSZF4B+Bdv/gKwI5afJuvYmknXLPVjQ5snpCAxJvSVtRnbkMio9RtImkrZuaJZsELLtJX0nzd8amB0RV5GNurpLw+sFtiArNGm2CicbsyQi5gIPpZPel7SwmeOBL0p6GphOy4cG/wFZF9ZUSdPT88bxvgu8LOkjBZMnkg1Stzg97pN+ExH3A7cAkyQ9QzZ2T9eIeBb4HnC/pKlkR0W9CraznOzE/4Hp0vD9ycZmeYrsXM9ladFdgUdi5UiiZiu46rNZHZP0CWDXiPheDcRyGTAuIv5a7Vis9vicjVkdi4i7aujqr2lONNYUH9mYmVnufM7GzMxy52RjZma5c7IxM7PcOdmYmVnunGzMzCx3TjZmZpY7JxszM8udk42ZmeXOycbMzHLnZGNmZrlzsjEzs9w52ZiZWe6cbMzMLHdONmZmljsnGzMzy52TjZmZ5c7JxszMcudkY2ZmuataspEUkm4qeN5O0hxJ9zSz3nBJP29mme6SvlJiHA+XFnGz7fSTNK0SbZmZrWuqeWTzLrCTpM7p+SHAG82tFBGPR8QZzSzWHSgp2UTE3qUsZ2ZmLVftbrTxwFHp8Wjg1oYZknaXNEnSU5IeljQwTd+/4ehH0gWSrpU0QdIMSQ1J6GJgW0lTJF0iqYukv0p6UtIzkkYWbGdRQbsTJI2V9E9JN0tSmrerpL9LekLSfZJ6FUx/WtLTwFfz3VVmZvWr2snmNuBYSZ2AwcCjBfP+CewTEcOA7wP/00Qb2wOHAbsD50tqD5wDvBwRQyPibGAJ8ImI2AU4APhJQyJpZBhwFjAI2AYYkdr7BTAqInYFrgV+mJa/Djg9Ioa06NWbma0n2lVz4xExVVI/sqOa8Y1mdwOulzQACKB9E83cGxHvA+9LegvYvMgyAv5H0r7Ah0DvtNy/Gy33WETMApA0BegHzAd2Av6c8lNb4E1J3YHuEfFgWvdG4IhmX7SZ2XqoqskmGQdcCuwP9CiY/gPggYj4REpIE5pY//2Cx8sp/pqOBzYFdo2IpZJmAp1KbEvA9IjYq3DBlGzMzKwE1e5Gg6xb6sKIeKbR9G6svGDgpDLbXAh0bdTWWynRHABsXUZbzwObStoLQFJ7STtGxHxgvqSPpuWOLzNGM7P1RtWTTUTMiohilzL/L/AjSU9R5hFYRMwFHpI0TdIlwM3AcEnPACeQnQ8qta0PgFHAj9OFAFOAhivYvgBcnrrcip0DMjMzQBFR7RjMzGwdV/UjGzMzW/c52ZiZWe6qnmwKytRc3Gj6uWW0cbWkQWuYP0HS8BbGt2+6GXSZpFGN5v1J0vzmSuy0NkkD0w2tDT/vSDqrBuLqJOmxdCPsdEkXVjsmgHRj8Fu1VG6oFmMCx1Uux7VS1ZMNWZmaF4BPN7rRsqRkI6ltRHwpIp7NJTp4jexquFuKzLsE+HxO222xiHg+3dA6FNgVeA+4q7pRAdml5Qemm2CHAodL2rO6IQEwBji82kE0MobaiwkcV7nG4LiA2kg2o4HLyD7UGy4vvhjonL6V39x4BUmLJP0kXR22V8ORi6S2ksakq9CekfT1Ruu1SfP/u9TgImJmREwluxm08by/kl1mXcsOIqum8Gq1A4nMovS0ffqp+hUq6cbc/1Q7jkK1GBM4rnI5rpWqelNnKlNzMPBfZMUzRwMPR8Q5kr6WvpkXsyHwaER8M7XTMH0o0DsidkrTuxes047sEuhpEfFD1h/HUlBzrtoktQWeAD4CXB4RjzazipmtA6p9ZHM0WZWAxcAdwMfTh1FzlqflG5sBbCPpF5IOB94pmHcF61mikdQBOAb4XbVjaRARy9OXiD7A7pJ2qnJIZtYKqp1sRgMHp/IxT5CVqzmwhPWWRMTyxhMjYh4whKy0zZeBqwtmPwwckI6m1hdHAE9GxOxqB9JYqsDwALXZn21mFVbNwdM2AvYBtoqIfhHRj6xM/+i0yNJUcbmcNnsCbSLiDuB7wC4Fs68hK/Z5u6RaqAnXGlYZtqHaJG3a0LWpbByjQyijmoOZ1a9qHtl8Avhbqtjc4PfAxyR1BK4Epha7QGANegMTUvmYm4DvFs6MiJ8CTwE3SirptUvaTdIs4NPAFZKmF8ybSNZFdZCkWZIOKyPWXEnakOzD/M5qx1KgF/CApKnAZODPEVH1y8Yl3QpMAgamv+MXHVNxjqs8jqtgmy5XY2Zmeav2ORszM1sPONmYmVnunGzMzCx3LUo2kraQdJuklyU9IWm8pO0qHVylSJqZrlSrdLvflfSSpOdr7OKAmqvHJKmvpAckPZvqop1Z7Zigpuu1Oa4yOK46iCkiyvohGyRsEvDlgmlDgH3KbWttfoB2ZSw7E+hZ4e0PAp4GOgL9gZeBtq25D9YQ275kl31Pq3YsBTH1AnZJj7uS1cMbVANxCeiSHrcHHgX2dFyOa12Nq1oxteTI5gBgaUT8pmFCRDwdERMBJJ0tabKkqQ0ZU1I/Sc9Juipl0vvTfRZIOiN9250q6bY0bRNJd6dpj0ganKZfIOlGSQ+RXb68qaQ70vYmSxqRluuRtjFd0tXkM4rmSOC2iHg/Il4BXgJ2z2E7ZYsarMcUEW9GxJPp8ULgObJL1asqMrVYr81xlcFxla5aMbUk2exEdrf/aiQdCgwg+9AdCuwqad80ewBZLawdgfnAp9L0c4BhETGY7K5/gAuBp9K0c4EbCjYzCDg4IhoKeP4sInZL7TVUDDgf+Efa1l3AVi14nc3pDbxe8HwWNfDhWQ8k9QOGkX2jqjplBVynAG+R3fvjuNbAcZWnFuOqRkyVvkDg0PTzFPAksD1ZkgF4JSKmpMdPAP3S46nAzZI+ByxL0z4K3AgQEX8DeiirOAAwLrJaapAV8fxl2mnjgI0kdSHrRroprX8vMK+ir9JaLP197gDOioh3mlu+NUSN1mtzXOVxXKWrRkwtSTbTycZIKUbAjyKNpRIRH4mIa9K8wkoBy1lZcfoo4HKycwyT1XwpmXcLHrch62ts2F7vgsPDvL0B9C143idNsyYoKz90B3BzRNRSZQOgduu1Oa7yOK7StWZMLUk2fwM6Sjq1YYKkwZL2Ae4DTk7fXpHUW9JmTTWkrGRM34h4APgO0A3oAkwEjk/L7A+83cS34PuB0wvaG5oePggcl6YdAWzcgtfZnHHAsZI6SupPdgT3WA7bWSdIEll9uuciKxtUE1Sj9docV3kcV+3HVHZByogISZ8A/k/Sd4AlZFd7nRURL0raAZiUfbawCPgc2ZFMMW2BmyR1Izsq+nlEzJd0AXCtshpa7wEnNrH+GcDlabl2ZEnmy2TnfG5VVsfsYbKB2SoqIqZLuh14lqz776tRpBJ1NSire7Q/0FNZXbfzC44wq2UE2aimz6RuT4BzI2J89UICsqvkrlc2tEUb4PaogXptOK5yOa4aj8m10czMLHeuIGBmZrlzsjEzs9w52ZiZWe5aLdlImqCshtiU9DO2wu33k3RcJdssYZu1Whut5uoxAUjqLmmspH8qqyixVw3ENLDgf3KKpHckneW4HNe6Gle1Ymq1CwQkTQC+FRGP59T+/qn9o/Nov8j2BpENubw7sCXwF2C7WrgiLV1mvGFELEr3tvwDODMiHqlyXNcDEyPiakkdgA3Sdf41IV2d8wawR0S8Wu14Gjiu8jiu0rVmTFXtRpPUTdKr6X4bJG0o6XVJ7SVtK+lPyqpKT5S0fVpmjKSfS3pY0gxJo1JzFwP7pEz9dUk7pm/3U5TVWBvQVBwtVMu10WquHlO6vH1fsnttiIgPainRJAcBL9fKB0EBx1Uex1W6VouptZPNzQWHbpdExAJgCrBfmn80cF9ELAWuBE6PiF2BbwG/KminF1lJm6PJkgxkNdYmpkoCPyO73+ayVJJhOFntskqq6dpoqr16TP2BOcB1kp6SdLWkDascU2PHkh2t1hrHVR7HVbpWi6m1k83xBaVlzk7Tfgt8Nj0+FvitsgoEewO/Sx+YV5AlmAZ3R8SHEfEssHkT25oEnKvsxtOtC+qprRdqsB5TO7KSRL+OiGFkZYfOqW5IK6VuvWOA31U7lkKOqzyOq3StHVMtXI02Djhc0iZkNdf+RhbX/ILENDQidihYp7DOWtHhAyLiFrIduRgYL+nACsddF7XRaqge0yxgVsER1liy5FMrjgCejIjZ1Q6kEcdVHsdVulaNqerJJp1bmEw2XMA96Rv5O8Arkj4N2QlvSUOaaWoh2aBcpHW2AWZExM+B3wODKxx6zdZGUw3WY4qIfwOvSxqYJh1EVuqnVoym9ro4wHGVy3GVrlVjau2r0XqRHWlAVlzz4DRvFNmh3P4R8fc0rT/w67ROe7KT8RdJGkOWlMam5RZFRJd01dV9QA9gDNkImp8HlgL/Bo6LiIoOKCbpPOBkstpoZ0XEHyvZfkspG2zuerLacw21jy6qblQrCqVeDXQAZgBfiIiqD/+Qzh29BmyTziPWBMdVHsdVumrE5NpoZmaWu6p3o5mZ2brPycbMzHLnZGNmZrlzsjEzs9w52ZiZWe6cbMzMLHdONmZmljsnGzMzy52TjZmZ5c7JxszMcudkY2ZmuXOyMTOz3DnZmJlZ7pxszMwsd80mG0kh6ScFz78l6YJm1vm4pEEViK9Y21c313ap25f0ZUknVCiuMWlcHjMza6SUI5v3gU9K6llGux8Hckk2EfGliGhuhMeSth8Rv4mIGyoSmJmZNamUZLMMuBL4euMZkvpJ+pukqZL+KmkrSXsDxwCXSJoiadtG64yR9GtJj0iaIWl/SddKei6Nwtmw3K8lPS5puqQLC6ZPkDQ8PV4k6YeSnk7tbV5s+5JOkTQ5LXeHpA3S+hdI+lZBuz+W9JikFyTtk6a3lXRJWn+qpP9K0yXpl5Kel/QXYLMy9ruZ2Xql1HM2lwPHS+rWaPovgOsjYjBwM/DziHgYGAecHRFDI+LlIu1tDOxFlsDGAT8DdgR2TkMHA5wXEcOBwcB+aajjxjYEHomIIcCDwClNbP/OiNgtLfcc8MUmXme7iNgdOAs4P037IrAgInYDdgNOSUNWfwIYSHYEdQKwdxNtmpmt90pKNhHxDnADcEajWXsBt6THNwIfLXG7f4hsPOpngNkR8UxEfAhMB/qlZT4j6UngKbJEVKxb7APgnvT4iYJ1G9tJ0kRJzwDHp/aKubNIW4cCJ0iaAjwK9AAGAPsCt0bE8oj4F/C3Nb1gM7P1Wbsylv0/4Engugps9/30+8OCxw3P26Ujh28Bu0XEvNS91qlIO0tT0gJYTtOvZwzw8Yh4WtJJwP7NxFXYloDTI+K+wgUlHdlEG2Zm1kjJlz5HxH+A21m1C+ph4Nj0+HhgYnq8EOi6FnFtBLwLLJC0OXBEmes33n5X4E1J7VOc5bgPOC2ti6TtJG1I1m332XROpxdwQJntmpmtN8q9z+YnQOFVaacDX5A0Ffg8cGaafhtwtqSnGl8gUIqIeJqs++yfZN10D5XZROPt/z+yLrCHUpvluBp4FnhS0jTgCrKjnruAF9O8G4BJZbZrZrbe0MpeKDMzs3y4goCZmeXOycbMzHJX9WQjqZ2kOZIubjT93DLaWGMJm8IbQVsQ376SnpS0rLAcjaShkialm06nSvpsS9rPi6TuksZK+me6YXavKsfTV9IDkp5N++zM5tdqHZIOTzfnviTpnGrH08BxlacW46rFmKBKcUVEVX/IrjR7CHiZdA4pTV9U4vptS1hmAjC8hfH1I7ux9AZgVMH07YAB6fGWwJtA92rvz4L4rge+lB53qHZsQC9gl/S4K/ACMKgG9lPb9L+3TdpPTzsux7WuxlTNuKp+ZAOMBi4DXiO7SZR0lNM5lZu5ufEKqUzNTyQ9DezVcOSSLkMeI2mapGckfb3Rem3S/P8uNbiImBkRU8nuASqc/kJEvJge/wt4C9i0vJeej1TpYV/gGoCI+CAi5lczpoh4MyKeTI8XklVy6F3NmJLdgZciYkZEfEB2JePIKscEjqtctRhXLcYEVYqrqslGUifgYOAPwK1kiYeIOAdYHFm5mWL3xWwIPBoRQyLiHwXThwK9I2KniNiZVW9AbUdWUufFiPhehV/H7mTfEIqV5qmG/sAc4Lp0+ffV6d6gmiCpHzCM7HL0ausNvF7wfBa1kQQdV3lqMa5ajAmqFFe1j2yOBh6IiMXAHcDHJbUtYb3lafnGZgDbSPqFpMOBdwrmXQFMi4gfrm3QhdINnTcCX4is5E4taAfsAvw6IoaR3SBbE/3FkrqQ/e3OiqwMkpmtB6qdbEYDB0uaSVaPrAdwYAnrLYmI5Y0nRsQ8YAjZOZovk92Q2eBh4IB0NFURkjYC7iUrGvpIpdqtgFnArIhoOHIYS5Z8qipVYbgDuDki7mxu+VbyBtC34HmfNK3aHFd5ajGuWowJqhRX1ZJN+qDeB9gqIvpFRD/gq6SuNGBpQ4mYMtrsCbSJiDuA77HqB+w1wHjgdknl1IRralsdyKoI3BARY9e2vUqKiH8Dr0samCYdRFbpoGokiexv8FxE/LSasTQyGRggqX/6mx5LVjW82hxXeWoxrlqMCaoVVxWviDgRuK3RtE3IzjV0BH5MdhL55iLrLmr0fAIwnOyo5klgSvo5onB+enwh2fmhNiXGuRvZkcK7wFxgepr+OWBpwbamAEOrtT+LxD0UeByYCtwNbFzleD4KRIqnYX8dWe39lGI7kuzquJfJjlKrHpPjWjfiqsWYqhWXy9WYmVnuqn3OxszM1gNONmZmljsnGzMzy12Lko2kLSTdJullSU9IGi9pu0oHVymSZqYr1Srd7ndTbaHnJR1W6fZbStK1kt5K4+/UDMdVulqMCRxXuRzXSmUnm3QJ613AhIjYNiJ2Bb4LbF7p4JqJY60vX17L7Q8iu2RwR+Bw4Fcl3pDaGsaQxVRrxuC4SjWG2osJHFe5xuC4gJYd2RwALI2I3zRMiIinI2IigKSzJU1WVgn5wjStn7LKw1cpq/h7v6TOad4ZyioBT5V0W5q2iaS707RHJA1O0y+QdKOkh4AbJW0q6Y60vcmSRqTleqRtTJd0NaC12UlNGEl26fb7EfEK8BJZzaGqi4gHgf9UO47GHFfpajEmcFzlclwrtSTZ7ER2t/9qJB0KDCD70B0K7Cpp3zR7AHB5ROwIzAc+laafAwyLiMFkd/1Ddi/MU2nauWQVlxsMAg6OiIYCnj+LiN1Sew0VA84H/pG2dRewVQteZ3Nqte6RmVnNqXRX1KHp56n0vAtZknkNeCUipqTpT5CV7ofsJr+bJd1NdvMhZDcAfgogIv6WjlQ2SvPGRVZLDbIinoOynj0ANkq1t/YFPpnWv1fSvMq9RDMzK1dLks10YFQT8wT8KCKuWGViVuX3/YJJy4HO6fFRZMnhY8B5knZuZvvvFjxuA+wZEUsaba+ZJiqiVusemZnVnJZ0o/0N6Cjp1IYJkgZL2ge4Dzg5HV0gqbekzZpqSFIboG9EPAB8B+hGdjQ0ETg+LbM/8HYUrxB8P3B6QXtD08MHgePStCOAjVvwOpszDjhWUkdJ/cmO4B7LYTtmZnWv7GQTWX2bT5BVa35Z0nTgR8C/I+J+4BZgkqRnyKoNd11Dc22Bm9KyTwE/j2yQrwvIzvdMBS4mq6NWzBnA8HQhwbOses5n3xTbJ8m68SoqIqYDt5MVuPwT8NUoUom6GiTdCkwCBkqaJemL1Y4JHFe9xwSOq1yOq2Cbro1mZmZ5cwUBMzPLnZONmZnlzsnGzMxy12rJRtKEVENsSvqp6OiWqUrBcZVss4Rt1mRtNFhRD+6ZtK8fr3Y8AJLOlDQtVXY4q5W3vVotKEmfTrF8KGl4a8bjuNaduGoxplqMq7WPbI6PiKHpp6l7dVqqH+ly59ag2q6N1uCAtK+r8s9eSNJOwClk1SWGAEdL+kgrhjCG1WtBTSO7WvHBVoyjsTE4rnKMofbiGkPtxQQ1FldVu9EkdZP0arrfBkkbSnpdUntJ20r6k7Kq0hMlbZ+WGSPp55IeljRDUkPSuhjYJ32T/7qkHSU9lp5PlTSgwuHXbG20GrUD8GhEvBcRy4C/k6o8tIZitaAi4rmIeL61YijGcZWnFuOqxZhSDDUVV2snm5sLutEuiYgFZGPR75fmHw3cFxFLgSuB01NV6W8BvypopxdZSZujyZIMZDXWJqZv8j8ju+fmsogYCgwnq11WSbVeGy2A+1OyPrXZpfM3jezLQA9JG5CNgd63mXXMbB3R2mX6j4+IxucPfgt8FniArFvqV8oqEOwN/E4rS890LFjn7oj4EHhWUlNDG0wiK3/TB7gzIl6s1IuoEx+NiDeUVXD4s6R/pm86VRERz0n6MVnVh3fJvmTUxE2wZpa/WrgabRxwuKRNgF3JyuG0AeYXnN8ZGhE7FKxTWGetaCG0iLgFOAZYDIyXdGCF467p2mgR8Ub6/RZZ5euqd/FFxDURsWtE7AvMA16odkxm1jqqnmwiYhEwmWy4gHsiYnmqg/aKpE9DNmCbpCHNNLWQgtI4krYBZkTEz4HfA4MrHHrN1kZL5766Njwmq8Rd9ZEC01EWkrYiO19zS3UjMrNWExGt8gNMAJ4n6z6ZAvylYN4osnMM+xVM609Wc+xpsvpj30/TxwCjCpZblH63Jzsqehr4Otk5nOlpW38CNsnhNZ0HvJxe1xGttS9LiGubtB+eTvvgvGrHlOKamP6WTwMHtfK2bwXeBJaSnV/7IlmNv1lkR8qzyc4XtvY+cVx1HlctxlSLcbk2mpmZ5a7q3WhmZrbuc7IxM7PcOdmYmVnunGzMzCx3TjZmZpY7JxszM8udk42ZmeXOycbMzHLnZGNmZrlzsjEzs9w52ZiZWe6cbMzMLHdONmZmljsnGzMzy52TjZmZ5c7JxszMcudkY2ZmuXOyMTOz3DnZmJlZ7v4/+D+jePGd8y0AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZsAAAGcCAYAAAAPnQ1AAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAABCqElEQVR4nO3debxVVf3/8deb2QEhBQ1BQwtRlEHFeZ6HDLKsHMqptMEhG8ypXw5l2VfLspwnnIc0jRTTTElTVFARAXNGxQzRQEFFGT6/P9a6eLicy70Xzz7nHng/H4/zuPusvfdan7PvOedz9tp7r62IwMzMrEjtah2AmZkt+5xszMyscE42ZmZWOCcbMzMrnJONmZkVzsnGzMwK52RjbYakQyX9q9ZxAEiaJGnHCtQzRdKunzyi5ZukkPS5JuaNlvStasdUa23p89ISTjY1kr+EPpA0W9I0SSMkrbyUdX2iD5ukvpLul/S+pH+31S9HSd/M8c3K22yUpK5FtBURG0bE6CLqXt7U05eipEskPStpgaRDy8xfV9Id+T34lqT/q0GYdcnJpra+EBErA5sAQ4GftmZlJZX4H94APAmsBpwC3CKpZwXqrRhJOwC/BA6IiK7ABsBNS1lXh0rGtiyQ1L7WMbQRTwHfA55oPENSJ+DvwH3Ap4E+wLVVja6OOdm0ARHxOnAXsJGkT+VfTtMlzcjTfRqWzXsxZ0p6CHgfuAbYDvhj3kv6o6TzJf2mtA1JIyX9oHHbktYjJbtTI+KDiLgVeBr4cktil/R5SU9KelfSa5JOK5nXN3d/HJbnzZD0HUmbSZogaaakP7ZwM20GjImIJ/M2+19EXBURs0q2y8K9u8a/pnMcR0l6Hnhe0oWSzmn0Wv4i6Yd5eoqkXSWtmfdAVy1ZbuP8q7ajpM9Kuk/S27nsOkndW/ialqi5uiVtkrf9LEl/knSTpF+UzP+JpDck/UfSt0q7ovKe9IV57/A9YKf8Wm/N772XJR1bUtcKkq7K/8Nnct1TS+afKOnFHMtkSfvm8g2Ai4Ct8vtzZi7vLOkcSa8q7aVeJGmFkvqOL4n98BZsrs9Keiy/D//S8P+SdKekYxpt1wkN8TUWEedHxD+AOWVmHwr8JyJ+GxHvRcSciJhQrh5Jp0v6Q57uKOk9SWeXbMs5JTFuKenh/Hl4SiXdt5K6Sbo8b4vXJf1CTfwwkHS2pH/ldT4n6Z+S3snvnaX6YVZREeFHDR7AFGDXPL0WMAn4OWnv4svAikBX4E/A7SXrjQZeBTYEOgAdc9m3SpbZHPgP0C4/70FKTGuUiWNf4JlGZX8E/pCntwVmLuF17AgMJP1wGQRMA76Y5/UFgvRl0wXYnfQhvh1YHegNvAnskJc/FPhXE+1sB3wAnA5sA3RuNL/xNlikrhzH34FVgRWA7YHXAOX5n8r1r1nm/3MfcERJXWcDF+XpzwG7AZ2BnsADwO/K/Z+X4j3SZN1AJ+AV4Pv5PfAl4CPgF3n+nsB/8/tkRdIv8AA+l+ePAN7J27JdXuZx4Ge57nWBl4A98vJnAf/M26kPMAGYWhLrV4A1c11fA94DejX1fwXOBUbm/0dX4K/Ar0pinwZsBKwEXF8ae5ntNBp4vWT5W4Fr87yvAo+WLDsYeBvo1My2/xdwaKOyK0g/7u4C3srtDmxi/Z2Bp/P01sCLDXHkeU/l6d45nr3zttstP++Z598GXJxf1+rAY8C3S7drXu9S4G5gxTzvBlIvRTvSZ2/bmn/n1TqA5fVB+hKaDcwkfWlcAKxQZrkhwIyS56OBMxotM5qSL9pc9gywW54+GhjVRBzfAB5pVHYmMGIpX9fvgHPzdN/8JdG7ZP7bwNdKnt8KHJenD6WJZJPn70X6UpqZt91vgfbltkHjunIcO5c8Fylpb5+fHwHc1+j/05BsvtUwL6/3WsN6ZWL8IvBkuXoq8J5ZWDcpWb5OTpa57F98nGyuIH955+efY/Fkc3XJ/C2AVxu1dxJwZZ5emHhKtsnUJcQ6HhjexP9CpGT02ZKyrYCXS2I/q2TeejSfbEqXH0BKvO1JX7QzgH553jnABS3Y1uWSzT3A3Pw+7AQcn7fLYomL9INmDunH44nAycBUYGXSD6bz8nInANc0Wvdu4BBgDeBDSr4XgAOA+0u266Ok7uRbS+MArgYuAfpU4r1XiYe70WrrixHRPSI+ExHfi4gPJK0o6WJJr0h6l/RrtnujXefXWlD3VcDX8/TXSb/IypkNrNKobBVgVktegKQtlE4umC7pHeA7pD2pUtNKpj8o87xFJ0ZExF0R8QXSr+HhpA9ba06MWLjdIn0ibyR9eAEOBK5rYr1bSd1AvUhf8guABwEkrSHpxtzF8S5pD6Lx61+MpLVzt9JsSbObWGZJda8JvJ5fx2KvL89/rYl55co+A6yZu3Jm5u6uk0lfeM3WJ+lgSeNL1t2IprdDT/KeVMnyf8vl5dp6pYl6mnotr5D29npExBzSl/HXlY5vHkDTn4XmfEBKmndFxEekxLUa6fjhIiLiA2AcsAPpPfNP4GHSnuQO+Tmk7f6VRtt9W6BXntcReKNk3sWkPZwGnyN9Fk7PMTX4CSmpP6Z0ZmVLuiIL5WTT9vwI6A9sERGrkN6okN44DRoP1V1u6O5rgeGSBpM+DLc30d4kYF0telbX4FzeEteTukPWiohupC4zLXmVTyYiFkTqV7+P9KUG6ZfyiiWLfbrcqo2e3wDsJ+kzpF/2tzbR3gzSr9qvkZLSjSVf8r/M9Q7M/6+v04LXHxGvRsTKDY8mFltS3W8AvSWVtrVWyfQbpO6ucvMWhlEy/Rppz6J7yaNrROzdXH15+11K2oNeLSK6AxNLYm283d8ifXFvWNJWt5Lt8EajeNcuE3tjjZefm9uB9MPrIGAX4P2IGNOC+sqZQPnPWlP+Seoy2xgYm5/vQermfiAv8xppz6Z0u68UEWfleR+SkmbDvFUiYsOSNp4BDgPuktS/oTAi/hsRR0TEmsC3gQvUxKnj1eJk0/Z0JX0QZ+YDiKe2YJ1ppD72hSJiKukNfg1wa/6ltZiIeI7U5XGqpC75wOkgmvjibSLe/0XEHEmbk76MK07ScEn7K51AodzWDsAjeZHxwJfynuHngG82V2ekkw3eAi4D7o6ImUtY/HrgYGC/PN2gK2nv8B1JvUldK5WypLrHAPOBoyV1kDSc9CXW4GbgMEkbSFoR+H/NtPUYMEvSCfkAdntJG0narKS+k/L2701KLA1WIn0JTweQdBgf/wiA9P7so3Q2FxGxgJSczpW0el6nt6Q9Sto6VNKAHHtLPgNfL1n+DOCWiJif2xtD2hv9Dc3s1UjqJKkLKVF2zJ+Jhu/Ja4EtlU4caQ8cR3r/PNNEdf8kvWcm572O0aQ98ZcjYnpJnV+QtEfe5l0k7SipT0S8QfqR8xtJq0hqp3TSyA6ljUTEDaS90HslfTa/jq/o4xOLZpD+PwuWuAUL5mTT9vyO1N/7FumL9G8tWOf3pF/oMySdV1J+FengfXPdBvuTTr2eQToQvF/Dh0HSdk1182TfA86QNIt0cPnmFsS7NGaQjqs8DzR0KZ0dEQ1dX+eS+umnkV53U11ijV0P7MqiCaSckUA/4L8R8VRJ+emks/neAe4E/tzCdluiybrzl9eXSEl1Jmmv5w7SL2Ei4i7gPOB+4AU+TsoflmsofzHvQzpG+DIfJ+FueZEzSMccXgbuBW4paWsy6Yt8DGn7DwQeKqn+PtKe8n8lNextnNAQV+4ivJe0R98Q++/yei/kv825hnQc6r+k4zTHNpp/dY6ruVOV7yH92NuadMzjA3LvQkQ8S9rOF5Hej8OBYY26r0o9TPosN+zFTCYdx2l4TkS8lus5mZSsXyP9qGj4bj6YdHxocm7zFlIX2yIi4irS/+g+SX1JZ28+mj+7I4HvR8RLzbz2QjWciWPLIEnbkz5cnwn/o5d5kh4lnSV3ZZl5G5C6tjpHxLwKtPVdYP+I2KHZhdsASQcDR0bEtrWOZXnlPZtllKSOpNNiL3OiWTZJ2kHSp3M32iGk7s+/lczfV+l6lk8Bvwb+urSJRlIvSdvkrpz+pGOLt1XidRQtd619j7SnYjXiZLMMyr9iZ5J2t39X02CsSP1JV7zPJH3575f7+Rt8m3Qd04uk4zvf/QRtdSKdCTWL1K31F9Lp+m1aPg40ndS911xXqRXI3WhmZlY479mYmVnhnGzMzKxwHv22jB49ekTfvn1rHYaZWV15/PHH34qIsiPGO9mU0bdvX8aNG1frMMzM6oqkJocWcjeamZkVzsnGzMwK52RjZmaF8zEbM6u6uXPnMnXqVObMKXdDTGvrunTpQp8+fejYsWOL16nrZCPpCtLggW9GxEZl5os0SOXepDtVHhoRi91b3Myqa+rUqXTt2pW+ffuy6F0SrK2LCN5++22mTp3KOuus0+L16r0bbQTpFrJN2Ys0Um8/4EjgwirEZGbNmDNnDquttpoTTR2SxGqrrdbqvdK6TjYR8QDwvyUsMpx069uIiEdId7xcbHhuM6s+J5r6tTT/u7ruRmuB3ix6u9ipueyN8ot/Mo9ccARdZzZ1H6VizO63L1t85UdVbdNsWSCJgw46iGuvTbe4mTdvHr169WKLLbbgjjvuaHK9cePGcfXVV3Peeec1uczMmTO5/vrr+d73vtdsHFtvvTUPP/xw619AI1OmTGGfffZh4sSJn7iuItT1nk0lSTpS0jhJ46ZPn978Ck1YMG9e1R5rffgiKz7zpwpuBbPlx0orrcTEiRP54IN0E9u///3v9O7du9n1hg4dusREAynZXHBBywbFrkSiqQfL+p7N6yx6b/I+uWwxEXEJ+X4XQ4cOXaqhsLf83qVLs9pSm/TLbWHeJ74Pltlya++99+bOO+9kv/3244YbbuCAAw7gwQcfBOCxxx7j+9//PnPmzGGFFVbgyiuvpH///owePZpzzjmHO+64g9NOO41XX32Vl156iVdffZXjjjuOY489lhNPPJEXX3yRIUOGsNtuu3HqqacyfPhwZsyYwdy5c/nFL37B8OHDAVh55ZWZPXs2o0eP5rTTTqNHjx5MnDiRTTfdlGuvvRZJPP744/zwhz9k9uzZ9OjRgxEjRtCrVy8ef/xxDj/8cAB23333mm3HlljWk81I0j3abwS2AN5pdL8PM6ux0/86icn/ebeidQ5YcxVO/cKGzS63//77c8YZZ7DPPvswYcIEDj/88IXJZv311+fBBx+kQ4cO3HvvvZx88snceuuti9Xx73//m/vvv59Zs2bRv39/vvvd73LWWWcxceJExo8fD6Quuttuu41VVlmFt956iy233JJhw4YtduzjySefZNKkSay55ppss802PPTQQ2yxxRYcc8wx/OUvf6Fnz57cdNNNnHLKKVxxxRUcdthh/PGPf2T77bfn+OOP/+QbrkB1nWwk3QDsCPSQNBU4FegIEBEXAaNIpz2/QDr1+bDaRGpmbdGgQYOYMmUKN9xwA3vvvfci89555x0OOeQQnn/+eSQxd+7csnV8/vOfp3PnznTu3JnVV1+dadOmLbZMRHDyySfzwAMP0K5dO15//XWmTZvGpz/96UWW23zzzenTpw8AQ4YMYcqUKXTv3p2JEyey2267ATB//nx69erFzJkzmTlzJttvvz0A3/jGN7jrrrs+8TYpSl0nm4g4oJn5ARxVpXDMbCm0ZA+kSMOGDePHP/4xo0eP5u23315Y/v/+3/9jp5124rbbbmPKlCnsuOOOZdfv3Lnzwun27dszr0zX9nXXXcf06dN5/PHH6dixI3379i176nC5uiKCDTfckDFjxiyy7MyZM1v5SmvLJwiY2XLt8MMP59RTT2XgwIGLlL/zzjsLTxgYMWJEq+rs2rUrs2bNWqSu1VdfnY4dO3L//ffzyitNDo68mP79+zN9+vSFyWbu3LlMmjSJ7t270717d/71r38BKaG1ZU42ZrZc69OnD8cee+xi5T/5yU846aST2HjjjcvurSzJaqutxjbbbMNGG23E8ccfz0EHHcS4ceMYOHAgV199Neuvv36L6+rUqRO33HILJ5xwAoMHD2bIkCELz2C78sorOeqooxgyZAipI6ftUlsPsBaGDh0a9XA/m0m/3JYF8+Yx8GeP1DoUs1Z55pln2GCDDWodhn0C5f6Hkh6PiKHllveejZmZFa6uTxAwCGDk+LKXDhWi6wod2an/6lVrz8yWDU429S6Cnl27VK256bM8JLyZtZ670czMrHBONmZmVjgnGzMzK5yTjZktlyTxox99fHuOc845h9NOO22J69x+++1Mnjy5kHi+9a1vNVt3S9u/6KKLuPrqqysS16GHHsott9zyietxsjGz5VLnzp3585//zFtvvdXidYpMNpdddhkDBgyoSPvf+c53OPjggysVWkU42ZjZcqlDhw4ceeSRnHvuuYvNmzJlCjvvvDODBg1il1124dVXX+Xhhx9m5MiRHH/88QwZMoQXX3xxkXUOPfRQvvvd77Lllluy7rrrMnr0aA4//HA22GADDj300IXLffe732Xo0KFsuOGGnHrqqQvLd9xxRxouJl955ZU55ZRTGDx4MFtuuSXTpk0r2/6ll17KZpttxuDBg/nyl7/M+++/D8Bpp53GOeecs7DeE044gc0335z11ltv4ajW8+fP5/jjj2ezzTZj0KBBXHzxxUAaNPToo4+mf//+7Lrrrrz55puV2d4VqcXMbGnddSL89+nK1vnpgbDXWc0udtRRRzFo0CB+8pOfLFJ+zDHHcMghh3DIIYdwxRVXcOyxx3L77bczbNgw9tlnH/bbb7+y9c2YMYMxY8YwcuRIhg0bxkMPPcRll13GZpttxvjx4xkyZAhnnnkmq666KvPnz2eXXXZhwoQJDBo0aJF63nvvPbbcckvOPPNMfvKTn3DppZfy05/+dLH2u3fvzhFHHAHAT3/6Uy6//HKOOeaYxeKaN28ejz32GKNGjeL000/n3nvv5fLLL6dbt26MHTuWDz/8kG222Ybdd9+dJ598kmeffZbJkyczbdo0BgwYsPCeOZ+E92zMbLm1yiqrcPDBBy92580xY8Zw4IEHAmno/obBLpvzhS98AUkMHDiQNdZYg4EDB9KuXTs23HBDpkyZAsDNN9/MJptswsYbb8ykSZPKdot16tSJffbZB4BNN9104bqNTZw4ke22246BAwdy3XXXMWnSpLLLfelLX1qsrnvuuYerr76aIUOGsMUWW/D222/z/PPP88ADD3DAAQfQvn171lxzTXbeeecWvfbmeM/GzGqrBXsgRTruuOPYZJNNOOywT367q4ZbBLRr126R2wW0a9eOefPm8fLLL3POOecwduxYPvWpT3HooYeWvdVAx44dF95YranbFkDqurv99tsZPHgwI0aMYPTo0UuMq7SuiOAPf/gDe+yxxyLLjho1qnUvuoW8Z2Nmy7VVV12Vr371q1x++eULy7beemtuvPFGIA3dv9122wGL3zqgtd59911WWmklunXrxrRp01p9s7PG7c+aNYtevXoxd+7cVt9iYI899uDCCy9ceFO45557jvfee4/tt9+em266ifnz5/PGG29w//33t6repjjZmNly70c/+tEiZ6X94Q9/4Morr2TQoEFcc801/P73vwfSbaTPPvtsNt5448VOEGiJwYMHs/HGG7P++utz4IEHss0227Rq/cbt//znP2eLLbZgm222adVtCyCdaj1gwAA22WQTNtpoI7797W8zb9489t13X/r168eAAQM4+OCD2WqrrVpVb1N8i4Ey6ukWA/PnzuW9rxez21vO9FlzGDakd9Xas2WTbzFQ/3yLATMza3OcbMzMrHBONmZmVjgnGzOrCR8vrl9L879zsjGzquvSpQtvv/22E04digjefvttunRp3U0bfVGnmVVdnz59mDp1KtOnT691KLYUunTpQp8+fVq1jpONmVVdx44dWWeddWodhlWRu9Hq3NxQrUMwM2uWk02dm7ug1hGYmTXPycbMzArnZGNmZoVzsjEzs8L5bLRlwBl3lL9hUhEG9e7mgTjNrNXqPtlI2hP4PdAeuCwizmo0f23gKqB7XubEiKjeMMkVcu7fn+P3/3h+kbIbO6WbID3zxuL31+ixcid6du28WPkn8crb7zN3ns9IMLPWq+tkI6k9cD6wGzAVGCtpZESU3mf1p8DNEXGhpAHAKKBv1YP9hH6w23r8YLf1Fi288gIeefltbjhiy6rEcMYdk5xszGyp1Psxm82BFyLipYj4CLgRGN5omQBWydPdgP9UMT4zM6PO92yA3sBrJc+nAls0WuY04B5JxwArAbuWq0jSkcCRAGuvvXbFAy3Sai//tSrtdPygG0S9/z4xs1qo92TTEgcAIyLiN5K2Aq6RtFFELNIfFBGXAJdAulNnDeJcanNX6FmVdqI9aO6HVWnLzJYt9f4z9XVgrZLnfXJZqW8CNwNExBigC9CjKtGZmRlQ/8lmLNBP0jqSOgH7AyMbLfMqsAuApA1IyWaZGWq2Z+f5tQ7BzKxZdZ1sImIecDRwN/AM6ayzSZLOkDQsL/Yj4AhJTwE3AIfGMnQTjdU7z6t1CGZmzar7Yzb5mplRjcp+VjI9Gdim2nGZmdnH6nrPxszM6oOTjZmZFa7uu9GsuoJg5PjGJ/wVo+sKHdmp/+pVacvMiuVkY63SsV07enbtUpW2ps+aU5V2zKx47kYzM7PCOdmYmVnhnGzMzKxwTjZmZlY4JxszMyuck42ZmRXOycbMzArnZGNmZoVzsjEzs8I52ZiZWeGcbMzMrHBONmZmVjgnGzMzK5yTjZmZFc7JxszMCudkY2ZmhXOyMTOzwjnZmJlZ4ZxszMyscE42ZmZWOCcbMzMrnJONmZkVzsnGzMwK52RjZmaFc7IxM7PCOdmYmVnh6j7ZSNpT0rOSXpB0YhPLfFXSZEmTJF1f7RjNzJZ3HWodwCchqT1wPrAbMBUYK2lkREwuWaYfcBKwTUTMkLR6baItRqePZtQ6BDOzZtX7ns3mwAsR8VJEfATcCAxvtMwRwPkRMQMgIt6scoyF6uxkY2Z1oN6TTW/gtZLnU3NZqfWA9SQ9JOkRSXtWLTozMwPqvButhToA/YAdgT7AA5IGRsTM0oUkHQkcCbD22mtXOUQzs2Vbve/ZvA6sVfK8Ty4rNRUYGRFzI+Jl4DlS8llERFwSEUMjYmjPnj0LC9jMbHlU73s2Y4F+ktYhJZn9gQMbLXM7cABwpaQepG61l6oZZNH6jvt5Vdr55QfwQLstgO2r0p6ZLTvqOtlExDxJRwN3A+2BKyJikqQzgHERMTLP213SZGA+cHxEvF27qJfS/b+Cf55VdtZKM55ZrOyjLj2Yu0Jl99DWWfAKxAJmV7RWM1seKCJqHUObM3To0Bg3blytw2iZ07oxabfqXDq04B8/T8nm63dVpb3ps+YwbEjj8z3MrK2S9HhEDC03r96P2ZiZWR1wsjEzs8LV9TEbq43VXv5rdRqa15l0boeZ1TsnG2sdtav4iQdNeus1Ro5vfCa7mRWp6wod2al/5Uf1crKpc//ptVutQyjMait1ga5dah2G2XJl+qw5hdTrYzZ17j9r7lHrEMzMmuVkY63y0YJaR2Bm9cjJxlplrpONmS0FJxszMyuck42ZmRXOycbMzArnU5+t1U58uDrt7NyzC1uuU522zKxYTjZW1nXPwvXPLVp2Y6f09+kyY2avvgKssWLl2n/pHdD8zmxZuSrNrIacbKysg/qnR6m+41KiufMLxbd/4sOkG0KY2TLBx2zMzKxwTjZmZlY4JxszMyuck42ZmRXOycZaZfUVah2BmdUjJxtrlUqe3mxmyw8nGzMzK5yTjZmZFc7JxszMCudkY2ZmhXOyMTOzwnlstDq3Qqf2zHz/o6q0NW/+AhZEVZoys2WMk02d23DNVWDlHtVp7IWOvPPB3Oq0ZWbLFHejmZlZ4ZxszMyscE42ZmZWOCcbMzMrXN0nG0l7SnpW0guSTlzCcl+WFJKGVjM+MzOr82QjqT1wPrAXMAA4QNKAMst1Bb4PPFrdCM3MDOo82QCbAy9ExEsR8RFwIzC8zHI/B34NzKlmcGZmltR7sukNvFbyfGouW0jSJsBaEXHnkiqSdKSkcZLGTZ8+vfKRmpktx+o92SyRpHbAb4EfNbdsRFwSEUMjYmjPnj2LD87MbDlS78nmdWCtkud9clmDrsBGwGhJU4AtgZE+ScDMrLrqPdmMBfpJWkdSJ2B/YGTDzIh4JyJ6RETfiOgLPAIMi4hxtQnXzGz5VNfJJiLmAUcDdwPPADdHxCRJZ0gaVtvozMysQd0PxBkRo4BRjcp+1sSyO1YjJjMzW1Rd79mYmVl9cLIxM7PCOdmYmVnhnGzMzKxwTjZmZlY4JxszMyuck42ZmRXOycbMzArnZGNmZoVzsjEzs8I52ZiZWeGcbMzMrHBONmZmVjgnGzMzK5yTjZmZFa7u72dj1bXSe6/Sd9zPC2/nlx/AA+22ALYvvC0zK56TjbXcujvy3gdzq/KmWWfBKxALmF2FtsyseE421nLr7cnEdkPpvmKnwpta8I+f89GCKLwdM1vUqKffYNiQ3hWv18dsrM2au0C1DsFsufO3SdMKqdfJxszMCudutHrXpRvMLuaXSDkrvvcGrLhu1dozs2WDk02967dbVZubP+XSqrZnZssGJxtr0864Y1KtQzCzCnCysZq77lm4/rlFy27MJ7w988asxZbvsXInenbtXIXIzJZd02d9yFuzPyo7r++Jdy5W9v1d+vGD3dZb6vacbKzmDuqfHqX6joOn34YbjtiyNkGZLacOuPQRppz1+YrX67PRzMyscE42ZmZWOCcbMzMrnJONmZkVzsnG2qw1usyvdQhmy509N1yjkHrrPtlI2lPSs5JekHRimfk/lDRZ0gRJ/5D0mVrEaa23RpcFtQ7BbLmz98BehdRb18lGUnvgfGAvYABwgKQBjRZ7EhgaEYOAW4D/q26UZmZW18kG2Bx4ISJeioiPgBuB4aULRMT9EfF+fvoI0KfKMZqZLffqPdn0Bl4reT41lzXlm8BdhUZkZmaLWW5GEJD0dWAosEMT848EjgRYe+21qxiZmdmyr96TzevAWiXP++SyRUjaFTgF2CEiPixXUURcAlwCMHToUN8isgkrdGrPzPfLj6dUSfPmL4AFPkHAbFlR78lmLNBP0jqkJLM/cGDpApI2Bi4G9oyIN6sf4rJlwzVXgZV7FN/QCx2Z9d7c4tsxs6qo62M2ETEPOBq4G3gGuDkiJkk6Q9KwvNjZwMrAnySNlzSyRuGamS236n3PhogYBYxqVPazkuldqx6UmZktoq73bMzMrD442ZiZWeGcbMzMrHBONmZmVjgnGzMzK5yTjZmZFc7JxszMCudkY2ZmhXOyMTOzwjnZmJlZ4ZxszMyscE42ZmZWOCcbMzMrnJONmZkVzsnGzMwK52RjZmaFc7IxM7PCOdmYmVnhnGzMzKxwTjZmZla4DrUOwOpMl24we1rx7cyfS6h98e2YWVU42Vjr9NutOu2Mu5L4aHp12jKzwrkbzczMCudkY2ZmhXOyMTOzwjnZmJlZ4XyCgLVZ7dqJ6bPm1DoMs+VK1xU6FlKvk421WV07d2DYkN61DsPMKsDdaGZmVjgnGzMzK5yTjZmZFa7uk42kPSU9K+kFSSeWmd9Z0k15/qOS+tYgTDOz5VpdJxtJ7YHzgb2AAcABkgY0WuybwIyI+BxwLvDr6kZpZmZ1nWyAzYEXIuKliPgIuBEY3miZ4cBVefoWYBdJqmKMZmbLvXpPNr2B10qeT81lZZeJiHnAO8BqjSuSdKSkcZLGTZ/uASBr7tMDYdV1ax2FmVWIr7PJIuIS4BKAoUOHRo3Dsb3OqnUEZlZB9b5n8zqwVsnzPrms7DKSOgDdgLerEp2ZmQH1n2zGAv0krSOpE7A/MLLRMiOBQ/L0fsB9EeE9FzOzKqrrbrSImCfpaOBuoD1wRURMknQGMC4iRgKXA9dIegH4HykhmZlZFdV1sgGIiFHAqEZlPyuZngN8pdpxmZnZx+q9G83MzOqAk42ZmRXOycbMzArnZGNmZoWTzwJenKTpwCtLuXoP4K0KhlO0eorXsRajnmKF+op3eYv1MxHRs9wMJ5sKkzQuIobWOo6Wqqd4HWsx6ilWqK94HevH3I1mZmaFc7IxM7PCOdlU3iW1DqCV6ilex1qMeooV6itex5r5mI2ZmRXOezZmZlY4JxszMyuck00FSdpT0rOSXpB0Yq3jKSVpLUn3S5osaZKk7+fy0yS9Lml8fuxd61gBJE2R9HSOaVwuW1XS3yU9n/9+qtZxAkjqX7L9xkt6V9JxbWXbSrpC0puSJpaUld2WSs7L7+EJkjZpA7GeLenfOZ7bJHXP5X0lfVCyfS9qA7E2+T+XdFLers9K2qOasS4h3ptKYp0iaXwur/y2jQg/KvAg3eLgRWBdoBPwFDCg1nGVxNcL2CRPdwWeAwYApwE/rnV8ZeKdAvRoVPZ/wIl5+kTg17WOs4n3wX+Bz7SVbQtsD2wCTGxuWwJ7A3cBArYEHm0Dse4OdMjTvy6JtW/pcm1ku5b9n+fP2lNAZ2Cd/F3RvtbxNpr/G+BnRW1b79lUzubACxHxUkR8BNwIDK9xTAtFxBsR8USengU8A/SubVStNhy4Kk9fBXyxdqE0aRfgxYhY2hEoKi4iHiDdy6lUU9tyOHB1JI8A3SX1qkqglI81Iu6JiHn56SOkO/LWXBPbtSnDgRsj4sOIeBl4gfSdUTVLileSgK8CNxTVvpNN5fQGXit5PpU2+mUuqS+wMfBoLjo6d1Fc0Va6poAA7pH0uKQjc9kaEfFGnv4vsEZtQlui/Vn0A9sWty00vS3b+vv4cNKeV4N1JD0p6Z+StqtVUI2U+5+39e26HTAtIp4vKavotnWyWc5IWhm4FTguIt4FLgQ+CwwB3iDtSrcF20bEJsBewFGSti+dGWlfv02dt59vTT4M+FMuaqvbdhFtcVuWI+kUYB5wXS56A1g7IjYGfghcL2mVWsWX1cX/vIwDWPRHUsW3rZNN5bwOrFXyvE8uazMkdSQlmusi4s8AETEtIuZHxALgUqq8a9+UiHg9/30TuI0U17SGLp38983aRVjWXsATETEN2u62zZralm3yfSzpUGAf4KCcHMldUm/n6cdJx0HWq1mQLPF/3ia3K4CkDsCXgJsayorYtk42lTMW6CdpnfwLd39gZI1jWij3yV4OPBMRvy0pL+2P3xeY2HjdapO0kqSuDdOkA8QTSdvzkLzYIcBfahNhkxb5ddgWt22JprblSODgfFbalsA7Jd1tNSFpT+AnwLCIeL+kvKek9nl6XaAf8FJtolwYU1P/85HA/pI6S1qHFOtj1Y6vCbsC/46IqQ0FhWzbap4Nsaw/SGfyPEf6FXBKreNpFNu2pK6SCcD4/NgbuAZ4OpePBHq1gVjXJZ258xQwqWFbAqsB/wCeB+4FVq11rCUxrwS8DXQrKWsT25aUAN8A5pKOFXyzqW1JOgvt/PwefhoY2gZifYF0vKPhfXtRXvbL+f0xHngC+EIbiLXJ/zlwSt6uzwJ7tYX3QS4fAXyn0bIV37YersbMzArnbjQzMyuck42ZmRXOycbMzArnZGNmZoVzsjEzs8I52ZhlkrpL+l7J8zUl3VKhuk+T9OM8fYakXStUby9Jd1Siribqn92KZe9tY0PyWBviZGP2se7AwmQTEf+JiP0q3UhE/Cwi7q1QdT8kXaneFlxDyfYzK+VkY/axs4DP5vt3nJ3v6TER0nApkm5XuvfLFElHS/phHqjwEUmr5uU+K+lveQDRByWt37gRSSMk7Zenp0g6XdITSvfvWT+Xr5QHcnwst9HUCOJfBv6W17lT0qA8/aSkn+XpMyQdkaePlzQ2DxR5eklMX89tjZd0ccPV4yXze0gaI+nzeW/qgbzsxJJBGkeSRlEwW4yTjdnHTiTdHmBIRBxfZv5GpDGkNgPOBN6PNFDhGODgvMwlwDERsSnwY+CCFrT7VqRBRy/M60C62vy+iNgc2Ak4Ow/ds1Ae9mRGRHyYix4EtpPUjTRg5Ta5fDvgAUm7k4Yd2Zw0UOSmkraXtAHwNWCbiBgCzAcOKmlnDeBO0r1O7gQOBO7Oyw4mXWVORMwAOktarQWv2ZYzHWodgFkduT/SvYBmSXoH+GsufxoYlEfU3hr4UxqKDkg3y2rOn/Pfx0nJDNJ4cMMajvMAXYC1SfchatALmF7y/EHgWOBlUnLYTdKKwDoR8Wzeu9kdeDIvvzIp+QwCNgXG5rhX4OOBOTuShrU5KiL+mcvGAlfkgV1vj4jxJTG8CaxJGrrHbCEnG7OW+7BkekHJ8wWkz1I7YGb+xb809c7n48+kgC9HxLNLWO8DUhJqMBYYShow8e9AD+AIUhJrqPNXEXFxaSWSjgGuioiTyrQxL6+/B/BPSDfhyrd8+DwwQtJvI+LqvHyXHJfZItyNZvaxWaRbZi+VSPcHelnSVyCNtC1p8FJWdzdwTB6tG0kbl1nmOdLtexva/4g0YOVXSF17D5K65R4oqfPwvAeGpN6SViftueyXp5G0qqTPNFRLumHZ+pJOyPM/Q7rR1qXAZaRbDTeMLP5p0i29zRbhZGOWRbp/x0P5oPfZS1nNQcA3JTWMWL20twb/OakLa4KkSfl543jfA16U9LmS4geBNyPigzzdJ/8lIu4BrgfGSHoauAXoGhGTgZ+S7ow6gbRX1KuknfmkA/8751PDdwSekvQk6VjP7/OimwKPxMe3cDZbyKM+m9UxSfsCm0bET9tALL8HRkbEP2odi7U9PmZjVsci4rY2dPbXRCcaa4r3bMzMrHA+ZmNmZoVzsjEzs8I52ZiZWeGcbMzMrHBONmZmVjgnGzMzK5yTjZmZFc7JxszMCudkY2ZmhXOyMTOzwjnZmJlZ4ZxszMyscE42ZmZWOCcbMzMrnJONmZkVzsnGzMwK52RjZmaFc7IxM7PC1SzZSApJ15Y87yBpuqQ7mllvqKTzmlmmu6TvtTCOh1sWcbP19JU0sRJ1mZkta2q5Z/MesJGkFfLz3YDXm1spIsZFxLHNLNYdaFGyiYitW7KcmZktvVp3o40CPp+nDwBuaJghaXNJYyQ9KelhSf1z+Y4Nez+STpN0haTRkl6S1JCEzgI+K2m8pLMlrSzpH5KekPS0pOEl7cwuqXe0pFsk/VvSdZKU520q6Z+SHpd0t6ReJeVPSXoKOKrYTWVmVr9qnWxuBPaX1AUYBDxaMu/fwHYRsTHwM+CXTdSxPrAHsDlwqqSOwInAixExJCKOB+YA+0bEJsBOwG8aEkkjGwPHAQOAdYFtcn1/APaLiE2BK4Az8/JXAsdExOClevVmZsuJDrVsPCImSOpL2qsZ1Wh2N+AqSf2AADo2Uc2dEfEh8KGkN4E1yiwj4JeStgcWAL3zcv9ttNxjETEVQNJ4oC8wE9gI+HvOT+2BNyR1B7pHxAN53WuAvZp90WZmy6GaJptsJHAOsCOwWkn5z4H7I2LfnJBGN7H+hyXT8yn/mg4CegKbRsRcSVOALi2sS8CkiNiqdMGcbMzMrAVq3Y0GqVvq9Ih4ulF5Nz4+YeDQVtY5C+jaqK43c6LZCfhMK+p6FugpaSsASR0lbRgRM4GZkrbNyx3UyhjNzJYbNU82ETE1Isqdyvx/wK8kPUkr98Ai4m3gIUkTJZ0NXAcMlfQ0cDDpeFBL6/oI2A/4dT4RYDzQcAbbYcD5ucut3DEgMzMDFBG1jsHMzJZxNd+zMTOzZZ+TjZmZFa7myaZkmJqzGpWf3Io6LpM0YAnzR0saupTxbZ8vBp0nab9G8/4maWZzQ+y0FZJ+IGlSPpZ1Q76+qU3KF+u+WQ9DADnWYtRTrFBf8dYi1ponG9IwNc8BX2l0oWWLko2k9hHxrYiYXEh08CrpbLjry8w7G/hGQe1WlKTewLHA0IjYiHS90P61jWqJRgB71jqIFhqBYy3CCOonVqiveEdQ5VjbQrI5APg96Uu94fTis4AV8nAz1zVeQdJsSb/JZ4dt1bDnIqm9pBH5l/vTkn7QaL12ef4vWhpcREyJiAmki0Ebz/sH6TTretGBtF07ACsC/6lxPE3KF8v+r9ZxtIRjLUY9xQr1FW8tYq3pRZ25G2dX4NukwTMPAB6OiBMlHR0RQ5pYdSXg0Yj4Ua6noXwI0Dv/cm984WUH0inQEyPiTJYzEfG6pHNISf0D4J6IuKfGYZnZcqLWezb7kEYJ+AC4FfiipPYtWG9+Xr6xl4B1Jf1B0p7AuyXzLmY5TTQAkj4FDAfWAdYEVpL09dpGZWbLi1onmwOAXfPwMY+ThqvZuQXrzYmI+Y0LI2IGMJg0tM13gMtKZj8M7NSWD4oXbFfg5YiYHhFzgT/z8cWpZmaFquXN01YBtgPWjoi+EdGXNEz/AXmRuXnE5dbU2QNoFxG3Aj8FNimZfTlpsM+b8zGL5c2rwJaSVswnYuwCPFPjmMxsOVHLPZt9gfvyiM0N/gJ8QVJn4BJgQrkTBJagNzA6Dx9zLXBS6cyI+C3wJHCNpBa9dkmbSZoKfAW4WNKkknkPAn8CdpE0VdIerYi1qiLiUeAW4AngadL//pKaBrUEkm4AxgD987b9Zq1jaopjLUY9xQr1FW8tYvVwNWZmVrhaH7MxM7PlgJONmZkVzsnGzMwKt1TJRtKnJd0o6UVJj0saJWm9SgdXKZKm5DPVKl3vSZJekPRsWz45AOpu3KYukh6T9FQey+30Wse0JPUUr2MtRj3FCjWKNyJa9SDdJGwM8J2SssHAdq2t65M8gA6tWHYK0KPC7Q8AngI6ky6UfBFoX81t0Mp4tyedCj6x1rG0IFYBK+fpjsCjwJa1jmtZiNexOtZaxbs0ezY7AXMj4qKShPVURDwIIOl4SWMlTWjIlpL6SnpG0qU5i94jaYU871hJk/PyN+ayVSXdnssekTQol58m6RpJD5FOX+4p6dbc3lhJ2+TlVsttTJJ0GcXcRXM4cGNEfBgRLwMvAJsX0E5FRH2N2xQRMTs/7Zgfbfa0yXqK17EWo55ihdrEuzTJZiPS1f6LkbQ70I/0pTsE2FTS9nl2P+D8iNgQmAl8OZefCGwcEYNIV/0DnA48mctOBq4uaWYAsGtENAzgeW5EbJbraxgx4FTgX7mt24C1l+J1Nqc38FrJ86m5zCpAaVDV8cCbwN8jXSfUZtVTvI61GPUUK1Q/3kqfILB7fjxJunhwfVKSgTRUyvg8/TjQN09PAK5TGqdrXi7bFrgGICLuA1ZTGnEAYGSksdQgDcHyx7zBRgKrSFqZ1GV0bV7/TmBGRV+lFS4i5kcaiLUPsLmkjWoc0hLVU7yOtRj1FCtUP96lSTaTgE2bmCfgVxExJD8+FxGX53mlIwXM5+MRpz8PnE86njBWzQ8l817JdDtSP2NDe71Ldg2L9jqwVsnzPrnMKigiZgL3Uyf3CamneB1rMeopVqhevEuTbO4DOks6sqFA0iBJ2wF3A4fnvQsk9Za0elMVKQ0Zs1ZE3A+cAHQDVgYeBA7Ky+wIvBUR75ap4h7gmJL6huTJB4ADc9lewKeW4nU2ZySwv6TOktYh7cE9VkA7y518LK57nl6BdIO9f9c0qCWop3gdazHqKVaoTbytHpAyIkLSvsDvJJ0AzCGd7XVcRDwvaQNgjNI9ZmYDXyftyZTTHrhWUjfSXtF5ETFT0mnAFZImAO8DhzSx/rHA+Xm5DqQk8x3SMZ8blMYxe5g0CGVFRcQkSTcDk0ndf0dFmZGo2wqlsZB2BHoojfV2asleZ1vTC7hK6XYT7YCbI6It33q7nuJ1rMWop1ihBvF6bDQzMyucRxAwM7PCOdmYmVnhnGzMzKxwVUs2kkYrjSE2Pj9uqXD9fSUdWMk6W9BmPY2Ntpak+/NoDZMkfb/WMTVFUv+S98l4Se9KOq7WcTWlnuJ1rMWop1ihNvFW7QQBSaOBH0fEuILq3zHXv08R9ZdpbwBwA2m0hDWBe4H12uoZaZJ6Ab0i4glJXUkX1n4xIibXOLQlymfLvA5sERGv1Dqe5tRTvI61GPUUK1Qv3pp2o0nqJumVfL0NklaS9JqkjpI+K+lvSqNKPyhp/bzMCEnnSXpY0kuS9svVnQVsl7P0DyRtqDSq6XilMdb6NRXHUqq3sdHeiIgn8vQs4BnqY3idXYAX6+FDm9VTvI61GPUUK1Qp3monm+tKdtvOjoh3gPHADnn+PsDdETEXuAQ4JiI2BX4MXFBSTy/SkDb7kJIMpDHWHswjCZxLut7m93k4hqGkscsqqW7HRpPUF9iYNNJrW7c/aQ+yXtRTvI61GPUUK1Qp3lZf1PkJHVSmG+0m4Guk4RL2By5QGoFga+BP+eJQSEP5N7g9IhYAkyWt0URbY4BTJPUB/hwRz1fqRdSzvG1vJV2EW25UhjZDUidgGHBSrWNpiXqK17EWo55iherG2xbORhsJ7ClpVdKYa/eR4ppZMubZkIjYoGSd0nHWyt4+ICKuJ23ED4BRknaucNx1NzaapI6kRHNdRPy51vG0wF7AExExrdaBtFA9xetYi1FPsUIV4615sskDZ44l3S7gjjwS6bvAy5K+AqBkcDNVzQK6NjyRtC7wUkScB/wFGFTh0OtqbDSlXcTLgWci4re1jqeFDqC+uiPqKV7HWox6ihWqGG+1z0brRdrTgDS45q553n7An4AdI+KfuWwd4MK8TkfSwfgzJI0gJaVb8nKzI2Ll/Kv9bmA1YASp2+0bwFzgv8CBEVHRm4dJOgU4nDQ22nERcVcl668kSduSBjh9GliQi0+OiFG1i6ppklYijWm3bj6216bVU7yOtRj1FCtUP16PjWZmZoWreTeamZkt+5xszMyscE42ZmZWOCcbMzMrnJONmZkVzsnGzMwK52RjZmaFc7IxM7PCOdmYmVnhnGzMzKxwTjZmZlY4JxszMyuck42ZmRXOycbMzArXbLKRFJJ+U/L8x5JOa2adL0oaUIH4ytV9WXN1t7R9Sd+RdHCF4hqR78tjZmaNtGTP5kPgS5J6tKLeLwKFJJuI+FZETK5E+xFxUURcXZHAzMysSS1JNvOAS4AfNJ4hqa+k+yRNkPQPSWtL2hoYBpwtabykzzZaZ4SkCyU9IuklSTtKukLSM/kunA3LXShpnKRJkk4vKR8taWieni3pTElP5frWKNe+pCMkjc3L3Sppxbz+aZJ+XFLvryU9Juk5Sdvl8vaSzs7rT5D07VwuSX+U9Kyke4HVW7HdzcyWKy09ZnM+cJCkbo3K/wBcFRGDgOuA8yLiYWAkcHxEDImIF8vU9ylgK1ICGwmcC2wIDJQ0JC9zSkQMBQYBO0gaVKaelYBHImIw8ABwRBPt/zkiNsvLPQN8s4nX2SEiNgeOA07NZd8E3omIzYDNgCPyLav3BfqT9qAOBrZuok4zs+Vei5JNRLwLXA0c22jWVsD1efoaYNsWtvvXSPejfhqYFhFPR8QCYBLQNy/zVUlPAE+SElG5brGPgDvy9OMl6za2kaQHJT0NHJTrK+fPZeraHThY0njgUWA1oB+wPXBDRMyPiP8A9y3pBZuZLc86tGLZ3wFPAFdWoN0P898FJdMNzzvkPYcfA5tFxIzcvdalTD1zc9ICmE/Tr2cE8MWIeErSocCOzcRVWpeAYyLi7tIFJe3dRB1mZtZIi099joj/ATezaBfUw8D+efog4ME8PQvo+gniWgV4D3hH0hrAXq1cv3H7XYE3JHXMcbbG3cB387pIWk/SSqRuu6/lYzq9gJ1aWa+Z2XKjtdfZ/AYoPSvtGOAwSROAbwDfz+U3AsdLerLxCQItERFPkbrP/k3qpnuolVU0bv//kbrAHsp1tsZlwGTgCUkTgYtJez23Ac/neVcDY1pZr5nZckMf90KZmZkVwyMImJlZ4ZxszMyscDVPNpI6SJou6axG5Se3oo4lDmFTeiHoUsS3vaQnJM0rHY5G0hBJY/JFpxMkfW1p6q8mSd0l3SLp3/ki2q1qHVNTJO2ZL5h9QdKJtY6nOfUUr2MthmNtRkTU9EE60+wh4EXyMaRcPruF67dvwTKjgaFLGV9f0oWlVwP7lZSvB/TL02sCbwDda709m3ktVwHfytOd2mq8QPv8flg3x/kUMKDWcS0L8TpWx1qrWGu+ZwMcAPweeJV0kSh5L2eFPNzMdY1XyMPU/EbSU8BWDXsu+TTkEZImSnpa0g8ardcuz/9FS4OLiCkRMYF0DVBp+XMR8Xye/g/wJtCzdS+9evLoD9sDlwNExEcRMbOmQTVtc+CFiHgpIj4inV04vMYxLUk9xetYi+FYm1HTZCOpC7Ar8FfgBlLiISJOBD6INNxMuetiVgIejYjBEfGvkvIhQO+I2CgiBrLoBagdSEPqPB8RP63w69ic9Auh3NA8bcU6wHTgynxK+GX5eqG2qDfwWsnzqbmsraqneB1rMRxrM2q9Z7MPcH9EfADcCnxRUvsWrDc/L9/YS8C6kv4gaU/g3ZJ5FwMTI+LMTxp0qXxB5zXAYZGG3GmrOgCbABdGxMaki2bbdL+ymS07ap1sDgB2lTSFNB7ZasDOLVhvTkTMb1wYETOAwaRjNN8hXZDZ4GFgp7w3VRGSVgHuJA0a+kil6i3IVGBqRDyan99CSj5t0evAWiXP++Sytqqe4nWsxXCszahZsslf1NsBa0dE34joCxxF7koD5jYMEdOKOnsA7SLiVuCnLPplejkwCrhZUmvGhGuqrU6kUQSujohbPml9RYuI/wKvSeqfi3YhjX7QFo0F+klaJ2/n/UkjebdV9RSvYy2GY23GJ/7S/QT2Be6LiNKBOP8C/J+kzqR76EyQ9EQTx23K6U06JtGQRE8qnRkRv80Hyq+RdFBLur0kbUZKKp8CviDp9IjYEPgq6YD7anlwT4BDI2J8C2OthWOA6/Ib7CXgsBrHU1ZEzJN0NGlcuvbAFRExqcZhName4nWsxXCszfNwNWZmVrhaH7MxM7PlgJONmZkVzsnGzMwKt1TJRtKnJd0o6UVJj0saJWm9SgdXKZKm5DPVKl3vSXlsoWcl7VHp+itJ0hWS3lS6J0+bV0/xOtZi1FOsUF/x1iLWVicbSSKdnTU6Ij4bEZuSzvpao9LBNRNHLc+kQ2ngz/2BDYE9gQtaeEFqrYwgxVkvRlA/8Y7AsRZhBPUTK9RXvCOocqxLs2ezEzA3Ii5qKIiIpyLiQQBJx0saqzQS8um5rK/SKMOXKo2SfI+kFfK8YyVNzsvfmMtWlXR7LntE0qBcfpqkayQ9RDp9uaekW3N7YyVtk5dbLbcxSdJlgD7JRmrCcODGiPgwIl4GXiCNOdQmRcQDwP9qHUdL1VO8jrUY9RQr1Fe8tYh1aZLNRqSr/RcjaXegH+lLdwiwqaTt8+x+wPn5GpWZwJdz+YnAxhExiHTVP8DpwJO57GTSiMsNBgC7RkTDAJ7nRsRmub6GEQNOBf6V27oNWHspXmdz6mksJDOzmqp0V9Tu+fFkfr4yKcm8CrxccsHj46Sh+wEmkC40vB24PZdtS05GEXFf3lNZJc8bmcdSgzSI54DUswfAKpJWJl1s+aW8/p2SZlTuJZqZWWstTbKZBOzXxDwBv4qIixcplPoCpSMFzAdWyNOfJyWHLwCnSBrYTPvvlUy3A7aMiDmN2mumioqop7GQzMxqamm60e4DOks6sqFA0iBJ25GGPzg8710gqbek1ZuqKA8rs1ZE3A+cAHQj7Q09CByUl9kReCsi3i1TxT2kIVga6huSJx8ADsxle5GGmqm0kcD+kjpLWoe0B/dYAe2YmdW9ViebSOPb7EsarflFSZOAXwH/jYh7gOuBMZKeJo0s3HUJ1bUHrs3LPgmcl2/odRrpeM8E4CzgkCbWPxYYmk8kmMyix3y2z7F9idSNV1F5LKGbSYNZ/g04qtxI1G2FpBuAMUB/SVMlfbPWMS1JPcXrWItRT7FCfcVbi1g9NpqZmRXOIwiYmVnhnGzMzKxwTjZmZla4qiUbSaPzGGLj86Oid7fMoxQcWMk6W9Bm3YyNBgvHiHs6b/9xbSCexcZnkvSVPPLDAklDaxlfY/UUr2MthmNdetXeszkoIobkR1PX6iytvuTTnatB9Tc2WoOd8vZvCx+KESw+PtNE0hmED1Q9muaNoH7iHYFjLcIIHOtSqWk3mqRukl7J19sgaSVJr0nqKOmzkv6mNKr0g5LWz8uMkHSepIclvSSpIWmdBWyXf7X/QNKGkh7LzydI6lfh8OtqbLS2qNz4TBHxTEQ8W6OQlqie4nWsxXCsS6/ayea6km60syPiHWA8sEOevw9wd0TMBS4BjsmjSv8YuKCknl6kIW32ISUZSGOsPZh/tZ9Luubm9xExBBhKGruskupxbLQA7skJ/MhmlzYzq5BqD9N/UEQ0PlZwE/A14H5St9QFeQSCrYE/6eOhZzqXrHN7RCwAJktq6tYGY0jD3/QB/hwRz1fqRdSxbSPidaVRHf4u6d/514+ZWaHawtloI4E9Ja0KbEoaDqcdMLPk+M6QiNigZJ3ScdbKDoQWEdcDw4APgFGSdq5w3HU3NlpEvJ7/vkkaDdvdfmZWFTVPNhExGxhLul3AHRExP4+D9rKkr0C6YZukwc1UNYuSoXEkrQu8FBHnAX8BBlU49LoaGy0fD+vaME0anbvN31HQzJYREVGVBzAaeJZ0jGY8cG/JvP1IxxN2KClbhzTm2FOk8cd+lstHAPuVLDc7/+1I2it6CvgB6RjOpNzW34BVC3hNpwAv5te1V7W25VLGum7eNk/l7XJKG4jpBuANYC7pmNc3SePuTSXtvU4jHcOr+fart3gdq2Nta7F6bDQzMytczbvRzMxs2edkY2ZmhXOyMTOzwjnZmJlZ4ZxszMyscE42ZmZWOCcbMzMrnJONmZkVzsnGzMwK52RjZmaFc7IxM7PCOdmYmVnhnGzMzKxwTjZmZlY4JxszMyuck42ZmRXOycbMzArnZGNmZoVzsjEzs8L9f53+ErSfBu1oAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 2 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Chi2=3.396391, p=0.065339 for all events secure, exploiting aggregates\n",
      "Chi2=3.396390, p=0.065339 for all 161 time moments secure\n"
     ]
    }
   ],
   "source": [
    "import sys\n",
    "from kmsurvival import main\n",
    "sys.argv[1:] = ['-i2', '--plot-curves']\n",
    "await main()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Complete Run: 5 logrank tests + survival tables\n",
    "To run the demo with three parties on localhost, for instance, we add `-M3` as command line option and run [kmsurival.py](kmsurvival.py) outside this notebook using a shell command. The plots are not shown this way, so instead we print the survival tables:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Using secure fixed-point numbers: SecFxp64:32\n",
      "Dataset: aml, with 3-party split, time 1 to 41 (stride 4) months\n",
      "2020-10-31 15:34:40,846 Logrank test on all events in the clear.\n",
      "Chi2=2.675902, p=0.101878 for all events in the clear\n",
      "2020-10-31 15:34:40,862 Start MPyC runtime v0.7\n",
      "2020-10-31 15:34:42,922 All 3 parties connected.\n",
      "2020-10-31 15:34:42,953 Logrank test on own events in the clear.\n",
      "Chi2=0.510517, p=0.474915 for own events in the clear\n",
      "          removed  observed  censored  entrance  at_risk\n",
      "event_at                                                \n",
      "0.0             0         0         0         4        4\n",
      "3.0             1         1         0         0        4\n",
      "5.0             1         1         0         0        3\n",
      "8.0             1         1         0         0        2\n",
      "12.0            1         1         0         0        1\n",
      "          removed  observed  censored  entrance  at_risk\n",
      "event_at                                                \n",
      "0.0             0         0         0         4        4\n",
      "2.0             1         1         0         0        4\n",
      "3.0             1         1         0         0        3\n",
      "7.0             1         1         0         0        2\n",
      "11.0            1         1         0         0        1\n",
      "2020-10-31 15:34:43,117 Logrank test on aggregated events in the clear.\n",
      "Chi2=2.685357, p=0.101275 for aggregated events in the clear\n",
      "          removed  observed  censored  entrance  at_risk\n",
      "event_at                                                \n",
      "0.0             0         0         0        11       11\n",
      "4.0             3         2         1         0       11\n",
      "8.0             4         3         1         0        8\n",
      "12.0            3         2         1         0        4\n",
      "44.0            1         0         1         0        1\n",
      "          removed  observed  censored  entrance  at_risk\n",
      "event_at                                                \n",
      "0.0             0         0         0        12       12\n",
      "4.0             6         5         1         0       12\n",
      "8.0             3         3         0         0        6\n",
      "12.0            3         3         0         0        3\n",
      "2020-10-31 15:34:43,242 Optimized secure logrank test on all individual events.\n",
      "2020-10-31 15:34:43,242 Interval 1 (time 1 to 4) # observed events = 7\n",
      "2020-10-31 15:34:43,258 Interval 2 (time 5 to 8) # observed events = 6\n",
      "2020-10-31 15:34:43,258 Interval 3 (time 9 to 12) # observed events = 5\n",
      "2020-10-31 15:34:43,258 Interval 4 (time 13 to 16) # observed events = 0\n",
      "2020-10-31 15:34:43,258 Interval 5 (time 17 to 20) # observed events = 0\n",
      "2020-10-31 15:34:43,258 Interval 6 (time 21 to 24) # observed events = 0\n",
      "2020-10-31 15:34:43,258 Interval 7 (time 25 to 28) # observed events = 0\n",
      "2020-10-31 15:34:43,258 Interval 8 (time 29 to 32) # observed events = 0\n",
      "2020-10-31 15:34:43,258 Interval 9 (time 33 to 36) # observed events = 0\n",
      "2020-10-31 15:34:43,258 Interval 10 (time 37 to 40) # observed events = 0\n",
      "2020-10-31 15:34:43,258 Interval 11 (time 41 to 41) # observed events = 0\n",
      "Chi2=2.675904, p=0.101877 for all events secure, exploiting aggregates\n",
      "2020-10-31 15:34:44,526 Secure logrank test for all 41 time moments.\n",
      "Chi2=2.675904, p=0.101877 for all 41 time moments secure\n",
      "2020-10-31 15:34:49,331 Stop MPyC runtime -- elapsed time: 0:00:08.469154\n"
     ]
    }
   ],
   "source": [
    "!python kmsurvival.py -M3 -i2 --print-tables --collapse"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To try out other runs of the demo for yourself, remember to consult MPyC's help message, using the `-H` option:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "usage: kmsurvival.py [-H] [-h] [-C ini] [-P addr] [-M m] [-I i] [-T t] [-B b]\n",
      "                     [--ssl] [-L l] [-K k] [--no-log] [--no-async]\n",
      "                     [--no-barrier] [--no-gmpy2] [--output-windows]\n",
      "                     [--output-file] [-f F]\n",
      "\n",
      "optional arguments:\n",
      "  -H, --HELP            show this help message for MPyC and exit\n",
      "  -h, --help            show kmsurvival.py help message (if any)\n",
      "\n",
      "MPyC configuration:\n",
      "  -C ini, --config ini  use ini file, defining all m parties\n",
      "  -P addr               use addr=host:port per party (repeat m times)\n",
      "  -M m                  use m local parties (and run all m, if i is not set)\n",
      "  -I i, --index i       set index of this local party to i, 0<=i<m\n",
      "  -T t, --threshold t   threshold t, 0<=t<m/2\n",
      "  -B b, --base-port b   use port number b+i for party i\n",
      "  --ssl                 enable SSL connections\n",
      "\n",
      "MPyC parameters:\n",
      "  -L l, --bit-length l  default bit length l for secure numbers\n",
      "  -K k, --sec-param k   security parameter k, leakage probability 2**-k\n",
      "  --no-log              disable logging messages\n",
      "  --no-async            disable asynchronous evaluation\n",
      "  --no-barrier          disable barriers\n",
      "  --no-gmpy2            disable use of gmpy2 package\n",
      "\n",
      "MPyC misc:\n",
      "  --output-windows      screen output for parties i>0 (only on Windows)\n",
      "  --output-file         append output for parties i>0 to party{m}_{i}.log\n",
      "  -f F                  consume IPython's -f argument F\n"
     ]
    }
   ],
   "source": [
    "!python kmsurvival.py -H"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
